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ABSTRACT 

The effectiveness of active control on asymmet- 
ric flows around circular cones is investigated com- 
putationally using cone spinning and rotatory oscil- 
lation around its axis. The investigation uses the 
time-accurate solution of the unsteady, compressible, 
full Navier-Stokes equations with the implicit, up- 
wind, flux-difference splitting, finite-volume scheme. 
The present solutions are obtained under the locally- 
conical-fiow assumption in order to understand the 
flow physics using very fine grids for reasonable flow 
resolution at low computational cost. For all the com- 
putational solutions, a grid of 241x81x2 points in 
the wrap-around, normal and axial directions, respec- 
tively, is used. The grid is spinning or oscillating 
rigidly with the cone according to its motion and 
the kinematical and dynamical boundary conditions 
are modified accordingly. The computational appli- 
cations include the effects of uniform spinning rates 
and periodic rotatory oscillations at different ampli- 
tudes and frequencies on the flow asymmetry. 

INTRODUCTION 

The problems of prediction, analysis and control 
of asvmmetric vortical flows around slender pointed 
bodies are of vital importance to the dynamic stability 
and controllability of missiles and fighter aircraft. 
The onset of flow asymmetry occurs when the relative 
incidence (ratio of angle of attack to nose semi-apex 
angle) of pointed forebodies exceeds certain critical 
values. In addition to the relative incidence as one 
of the influential parameters for the onset of flow 
asymmetry, the freestream Mach number. Reynolds 
number and shape of the body-cross sectional area 
are also important parameters. 
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Substantial research efforts have recently been de- 
voted for eliminating or alleviating flow asymmetry 
and its corresponding side force. In the experimental 
area, several passive-control methods * * ^ and active- 
control methods 4-8 have been investigated. Compu- 
tational simulations have also been used to investigate 
the effectiveness of several passive-control methods 9 
and active-control methods 7 ’ l0, 11 . Various methods 
of passive control were demonstrated in the above ref- 
erences which include the use of vertical fins along 
the leeward plane of geometric symmetry, thin and 
thick side strakes with different orientations, and ro- 
tatable forebody tips which have variable cross sec- 
tion (from a circular shape at its base to an elliptic 
shape at its tip). It was shown by Kandil et al. 9 that 
side-strake control is more practical than vertical-fin 
control since the former was more effective over a 
wide range of angle of attack than the latter. More- 
over, side-strake control provided an additional lifting 
force. However, the effectiveness of side-strake con- 
trol terminated at very high angles of attack for the 
considered strake geometry and flow conditions. 

Various active-control methods have been used 
which include forebody blowing and movable fore- 
body strakes. The forebody blowing methods include 
forward blowing, normal blowing, aft blowing and 
tangential blowing. The main purposes of forebody 
blowing are to control flow separation on the fore- 
body and to create yawing forces and moments which 
can be utilized in controlling the body. 

In Ref. 12, the authors investigated the effective- 
ness of two methods of active control which include 
flow injection and surface heating. The study of 
flow-injection control covered normal and tangential 
injection. Moreover, a hybrid method of flow con- 
trol which combined surface heating and normal in- 
jection methods was also investigated. These active 
control methods were directed at either rendering the 
asymmetric vortical flow symmetric or rendering the 
surface-pressure distribution symmetric. 

Active control of asymmetric flows around slender 
pointed bodies using body spinning about its axis has 
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recently been investigated experimentally. Kruse 13 
investigated the effects of several spinning rates on 
the side force of a 10° half apex-angle cone. The re- 
sults were presented for spin tests at 58-60° angle of 
attack, 0.6 Mach number and 1 x 10 6 Reynolds num- 
ber (based on diameter). For the 58° angle of attack, 
the side force versus the roll angle was shown for 
four spin rates which varied from 60 rpm (revolu- 
tion per minute) to 400 rpm. For all values of spin 
rates, the side force changed direction in an irregular 
manner within each revolution of the cone and it was 
repeatable from one revolution to the next. More- 
over, the side-force variation within each revolution 
showed roughly three cycles, and there was an evi- 
dence that the amplitude of the side force decreases 
with the increase of spin rate. 

Fidler 14 used spinning of the nose, nose tip and 
a band of the body surface as active-control methods 
for alleviating asymmetric vortex effects on a tangent- 
ogive configuration. By rotating the nose, nose tip 
and a band of the body just aft of the nose, the wake 
pattern and the associated side forces and moments 
were cyclically altered. For the nose and nose-tip 
rotations, the peak-to-peak variations of the side force 
were decreased as the spin rate was increased. The 
results also showed that the average side force was 
constant throughout the spin range. However, by 
using the nose tip with three axial grit strips, the mean 
side force was brought to zero. 

Experimental and computational studies on the 
effects of rotation and rotatory oscillations on the 
vortex shedding behind circular cylinders have re- 
cently been conducted by several researchers' ? '* 7 . 
Coutanceau and Menard ' 4 reviewed earlier work and 
conducted experimental investigation on a circular 
cylinder undergoing steady rotatory and rectilinear 
motion. They concluded that if the ratio of rota- 
tional velocity to rectilinear velocity is greater than 2, 
then the Karman vortex street disappears. Taneda 16 
showed that if the cylinder is forced to undergo a har- 
monic rotatory oscillations at large values of ampli- 
tude and frequency then vortex shedding is eliminated 
and a symmetric flow can be generated. Chen 17 et 
al. integrated the velocity/vorticity formulation using 
an explicit finite-difference/pseudo-spectral technique 
and the Biot-Savart law to study the temporal devel- 
opment of two-dimensional incompressible, viscous 
flow around a circular cylinder undergoing steady ro- 
tatory' and rectilinear motion. Their computational 
results showed that rotation does not suppress vortex 
shedding for large values of the ratio of rotational ve- 
locity to rectilinear velocity. This conclusion is not 


in agreement with that of Cautanceau and Menard. 

In the present paper, we investigate the effec- 
tiveness of spinning and rotatory oscillation as ac- 
tive control methods to eliminate or alleviate the side 
forces due to vortex asymmetry for a 5°-semi-apex 
angle, circular cone. The investigation uses the time- 
accurate solution of the unsteady, compressible, full 
Navier-Stokes (NS) equations. The locally, conical- 
flow assumption is used to obtain all the present so- 
lutions since it provides excellent flow physics at 
substantially lower computational cost in compari- 
son with that required for the corresponding three- 
dimensional flow cases. 

FORMULATION 


Governing Equations: 

The vector form of the governing equations is de- 
veloped in terms of an inertial frame of reference, 
and hence there are no source terms on the right-hand 
side of the equations. Hence, the components of the 
flow-field vector [p,pV,pe] < are absolute quantities. 
This is unlike the earlier development of the gov- 
erning equations by the principal author of this paper 
(Ref. 18), where the equations are developed in terms 
of a non-inertial frame of reference (translating and 
rotation frame of reference) and source terms appear 
on the right-hand side of the equations. 

The conservative form of the dimensionless, 
unsteady, compressible, full NS equations in 
terms of time-dependent, body-conformed coordi- 
nates £*,£ 2 and is given by 

= 0: m = 1 - 3, s = 1 - 3 

* ( 1 ) 

= £ m (xi, t) ( 2 ) 

Q = = j[p,pui,pu- 2 ,pm,pe]' (3) 


dQ dEm 
~dt + 

where 


E m = inviscid flux in direction 

d£ m 1 ' 


pU m ,pttlUm + d\ C'p,PU2Um 


+ d-2t m p,puzUm + dsZ m p,(pe + p)U m 

~ ~dT 1 


(4) 
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(E v ) g = viscous and heat-conduction flux in f 

direction 

= j[0,dk^T kh d k CTk2^kC T kZ^ 

dkS s {unTkn - 9Jb)] , i fc = l-3, n = l-3 (5) 


df m 

Um = BkCuk + -jf ( 6 > 

and i s the grid speed. The three momentum 
elements of Eq. (5) are given by 






; j = i - 3 


(7) 


The last element of Eq. (5) is given by 

dkf s { u p T kp ~ 9k) = (dkt s dpC 


1 




,a(« 2 ) 


;p = 1 - 3 


( 8 ) 


(•)-i)P r *"' ae 1 J 

The reference parameters for the dimensionless form 
of the equations are L,a. x .,L/a <x ..,p x ai 'd for 
the length, velocity, time, density and molecular vis- 
cosity, respectively. The Reynolds number is defined 
as Rf - p. x .Y x .LIn oo, where L is the cone length. 
The pressure, p, is related to the total energy per unit 
mass. e. and density, p, by the gas equation 

]> — (* — 1 )p^f — —Untln^j . (9) 

The viscosity, //, is calculated from the Sutherland 

law ‘- + c\ 

+ ' C = 0.4317. (10) 


- ***(? 


t + c r 


and the Prandtl number Pr = 0.72. 

In Eqs. (1)-(10). the indicia! notation is used for 
convenience. The subscripts j. k and n are summation 
indices, the superscript or subscript s is a summation 
index and the superscript or subscript m is a free 
index. The range of j. k. n. s and m is 1-3, and 
d t = d/dxk. In Eqs. (1)-(I0), u„ is the Cartesian 
velocity component, U m the contravariant velocity 
component. t „ the Cartesian component of the shear 
stress tensor, q * the Cartesian component of heat flux 
vector, a the local speed of sound and M x the free- 
stream Mach number. 


Boundary and Initial Conditions 
and Grid Motion: 

Boundary conditions are explicitly implemented. 
They include inflow-outflow conditions, solid-bound- 
ary conditions and plane of geometric symmetry con- 
ditions. At the plane of geometric symmetry, periodic 
conditions are enforced. Since we are dealing with 
a supersonic flow, at the far-field inflow boundaries, 
freestream conditions are specified and all the five 
flow variables are extrapolated from the exterior val- 
ues according to the Riemann-invariant characteristic 
conditions. The conical shock is captured as part of 
the solution. At the far-field outflow boundaries, first- 
order extrapolation from the interior point is used. 

Since the cone is undergoing spinning or rotatory- 
oscillation motion at an angular velocity of ue x 
around the x-axis, where e x is a unit vector along 
the x-axis and u> is either uniform (for spinning mo- 
tion) or time dependent (for rotatory oscillation), the 
grid is moved with the same angular velocity as that 
of the body. The grid speed, and the metric co- 
efficient, are computed at each time step of the 
computational scheme. Consequently, the kinemati- 
cal boundary conditions at the inflow-outflow bound- 
aries and at the cone surface are expressed in terms of 
the relative velocities. For the dynamical boundary 
condition, §£ at the cone surface is no longer equal 
to zero. This condition for the present rotating cone 
is modified as 


dp 

dn 


2 

= - pa c ■ n = pr c u; — p 


df 

dt 


-rrr ) /r. 


cone 


(ii) 


where r c is the cone radius and a c is the acceleration 
of the cone surface. It should be noted that the 
other tangential acceleration does not contribute to 
the right-hand side of Eq. (13) since its scalar product 
with ft is zero. Finally, the boundary condition for the 
temperature is obtained from the adiabatic boundary 
condition, | cone — 0* 

The initial conditions correspond to the asymmet- 
rir flow solution without the cone rotation. 


COMPUTATIONAL SCHEME 

The implicit, upwind, flux-difference splitting, 
finite-volume scheme is used to solve the unsteady, 
compressible, full Navier-Stokes equations. The 
scheme uses the flux-difference splitting of Roe. The 
smooth flux limiter is used to eliminate oscillations 
at locations of large flow gradients. The viscous- 
and heat-flux terms are linearized in time and the 
cross-derivative terms are eliminated in the implicit 
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operator and retained in the explicit terms. The vis- 
cous terms are differenced using second-order accu- 
rate central differencing. The resulting difference 
equation is approximately factored and is solved in 
three sweeps in the and £ 3 directions. The 
computational scheme is coded in the computer pro- 
gram “FTNS3D”. 

For the locally-conical flow solutions, an axial 
station of x\ — 1.0 is selected and the components of 
the flowfield vector are forced to be equal between 
this axial station and another axial station in close 
proximity. This ensures that the flow variables are 
locally independent of the axial direction at x\ = 1.0. 

The method of solution consists of two steps. In 
the first step, the problem is solved for the asymmetric 
flow without rotating the cone. This solution repre- 
sents the initial conditions for the second step. In the 
second step, the cone spinning rate or rotatory oscil- 
lation is specified and the NS equations are solved 
accurately in time. At each time step, the body and 
the grid are rotated through an angle corresponding 
to the cone rotational velocity. The metric coeffi- 
cients and the grid speed are computed and the Roe 
flux-difference splitting scheme is used to obtain the 
solution. The computations proceed until periodic 
responses are obtained. 

COMPUTATIONAL APPLICATIONS 

A 5° semi-apex angle circular cone is used for the 
present applications. The freestream Mach number 
is 1.8 and the Reynolds number based on the cone 
length is 1 0\ A grid of 241x81x2 grid points in 
the wrap-around, normal and axial directions, respec- 
tively. is used for all the computational applications. 
The minimum grid size, Af : \ at the cone surface is 
of order 10 -4 . 

Initial Conditions, Cone at o = 20°, = 0: 

The cone is set at an angle of attack, a, of 20° and 
the spinning rate, is set equal to zero. Figure 1 
shows the total-pressure-loss (TPL) contours, surface- 
pressure (SP) coefficient versus angle 9 (measured 
from the windward plane of geometric symmetry) 
and streamlines for the locally-conical flow solution 
after 13.000 time steps. The solution is obtained 
using time-accurate stepping with At = 0.001. The 
solution shows that the flow is steady and asymmetric. 
This solution serves as initial conditions for the next 
applications with different spinning rates. 


Uniform Spinning Rate: 

There are several basic ideas behind the use of 
spinning to alleviate or possibly eliminate the asym- 
metry of vortices and their subsequent result of pro- 
ducing a side force on the cone. One of the ideas 
can be explained by considering the asymmetric so- 
lution of Fig. 1, which represents the initial condi- 
tion. By spinning the cone in the counter-clockwise 
direction, the speed of boundary-layer flow on the 
right-hand side of the cone is enhanced for resist- 
ing flow separation while the speed of the boundary- 
layer flow on the left-hand side of the cone is re- 
tarded for producing early flow separation. More- 
over, the spinning motion is adding either positive 
or negative vorticity in the flowfield. Hence, by se- 
lecting the appropriate spinning rate, the asymmetric 
vortices might be rendered symmetric. This is the 
first effect of spinning. The second effect of spin- 
ning is the increase of the pressure gradient normal 
to the body (f£|cone = P^ T con^- For small values 
of ijj, the effect of the pressure gradient will not be 
pronounced. However, for large values of cj, the ef- 
fect of the pressure gradient will be pronounced, the 
other idea behind using spinning as an active con- 
trol method is based on the experimental data where 
Kruse 13 and Fidler 14 found that spinning produced 
oscillatoiy side-force response. If the mean of the 
side force is zero then the average side force will be 
zero. In the next applications, we present the effects 
of constant spinning at tangential velocities, of 
±0.06 and ±0.2 which correspond to spinning rates, 
ft, of ±2,292 rpm and ±7,639 rpm for a cone of one 
meter long, respectively. 

1. Uniform Spinning at ±0.06: 

Figures 2 and 3 show the results for uniform 

® . . d £ 2 

counter-clockwise (CCW) spinning at Aft- = 0.06, 

and Figures 4 and 5 show the results for uniform 
clockwise (CW) spinning at 5^- = —0.06. All spin- 
ning cases start at the time step n = 13,001 and the 
solutions are obtained with At = 0.001. With the 
present surface speed of 0.06, the cone rotates one 
revolution in 9.163 dimensionless time, which cor- 
responds to 9,163 time steps. Figure 2 shows the 
side-force and lift coefficients (Cy and Cl) versus the 
number of time steps. It is observed that the force 
coefficients reach a periodic response very quickly 
and the period of oscillation is 9.163 dimensionless 
time, which is equal to the time required to rotate 
the cone one revolution. The Cy curve oscillates be- 
tween -0.00069 and -0.00053 with a mean value of 
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-0.00061 . The magnitude of the mean Cy is less than 
the magnitude of Cy without spinning which is equal 
to 0.00065. Thus, the CCW spinning reduces the 
magnitude of the side force on the average. The Cl- 
coefficient curve shows small-amplitude periodic re- 
sponse. In Fig. 3, snapshots of the TPL contours, SP 
coefficient and streamlines are shown at five instants 
covering one cycle of periodic side-force response. 
They are marked by the numbers 1, 2, 3, 4 and 5 
on Figs. 2 and 3. The snapshots show that the right 
and left vortices heights, lateral positions, strength- 
ens and separation points oscillate slightly. Conse- 
quently, the corresponding surface pressures oscillate 
slightly too. The CCW spinning delays flow separa- 
tion on the right side and expedites flow separation 
on the left side. 

With the CW spinning of -0.06, the Cy and Cl 
curves of Fig. 4 show that their periodic response is 
also reached very quickly. The Cy curve oscillates 
between -0.00077 and -0.00055 with a mean value of 
-0.00066. The magnitude of the mean Cy is slightly 
higher than the magnitude of Cy without spinning. 
Thus, the CW spinning does not reduce the mean 
value of side force. Figure 5 shows snapshots of the 
TPL contours, SP coefficients and streamlines at four 
instants (marked by 1, 2, 3 and 4 in Figs. 3 and 4) dur- 
ing one cycle of periodic response. The CW spinning 
increases flow separation on the right side and delays 
flow separation on the left side. Comparisons of the 
snapshots of Fig. 5 with the corresponding snapshots 
of Fig. 3, show that the vortex on the right-hand side 
of Fig. 5 (snapshot 1) moves more to the right while 
the vortex on the right-hand side of Fig. 3 (snapshot 
1) moves more to the left. Similar motions are ob- 
served for the vortex on the left side of Figs. 5 and 3 
(snapshot 1 ). Hence, the side force at point 1 of the 
CW spinning will be of higher magnitude than the 
side force at point 1 of the CCW spinning. 

2. Uniform Spinning at ±0.2: 

Next, the uniform spinning is increased to 0.2. 
The results of the CCW spinning are shown in Figs. 6 
and 7 and the results of the CW spinning are shown 
in Fig. 8. With the surface speed of 0.2, the cone 
rotates one revolution in 2.749 dimensionless time, 
which corresponds to 2,749 time steps. Figure 6 
shows that the Cy and Cl curves reach a periodic 
response very quickly and the period of oscillation is 
2.749 dimensionless time. The Cy curve oscillates 
between -0.00089 and -0.000050 with a mean value 
of -0.00047. The magnitude of the mean Cy is 
substantially lower than the magnitude of Cy without 


spinning. Thus, the high CCW spinning does reduce 
the mean value of side force. It should be noticed 
that the amplitude of oscillation of the Cl curve is 
higher than that of Fig. 2. Three snapshots of the 
TPL contours, SP coefficient and streamlines during 
a half-cycle of the periodic Cy curve (marked by 1, 

2 and 3 in Figs. 6 and 7) are shown in Fig. 7. It is 
noticed that the CCW spinning substantially increases 
the flow separation on the left side and delays the 
flow separation on the right side. Also, it is noticed 
that the right-hand side vortex moves more in the 
downward and leftward directions than that of Fig. 3. 
The left-hand side vortex moves more in the leftward 
direction than that of Fig. 3. 

Figure 8 shows the Cy and Cl periodic responses 
for the CW spinning case at -0.2. The Cy curve 
oscillates between —0.00102 and —0.00021 with a 
mean value of -0.00615. The magnitude of the 
mean Cy is lower than the magnitude of Cy without 
spinning but it is substantially higher than the mean 
value of the CCW spinning of Fig. 6. A snapshot of 
the TPL contours, SP coefficient and streamlines is 
shown in Fig. 8. The CW spinning is observed to 
increase the flow separation on the right side and the 
left-hand side vortex moves more to the right. 

3. Uniform Spinning at +0.6: 

The uniform spinning is increased to 0.6 and the 
results of the CCW spinning are shown in Fig. 9. 
With the speed of 0.6, the cone rotates one revo- 
lution in 0.916 dimensionless time. The Cy curve 
shows the periodic response which oscillates between 
-0.005 and +0.0038 with a mean value of -0.0006. 
With this high value of CCW spinning, the side-force 
coefficient is oscillating between positive and nega- 
tive values and the vortices on the left and right sides 
are changing heights periodically. It is noticed that 
the boundary layer at certain instants will become a 
free-shear-layer band around the body (e.g. Fig. 9- 
snapshot 1). Although the mean value of side force 
is not zero, this test tells that there is a certain CCW 
spinning value at which the mean side force will be 
zero. 

Rotatory Oscillations: 

In this section, we investigate the effect of pe- 
riodic rotatory oscillation of the cone on the flow 
asymmetry and the side force. The form of the sur- 
face speed is given by 

d( 2 2 z _ 

— = V s cos — 1 (12) 

at t 
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where V, is the surface-speed amplitude and r is the 
period of oscillation. Substituting Eq. (12) into the 
relation jjf = ^ £/r c and integrating the result one 
obtains the corresponding equation for the angular 
motion, 6, which is 

Oir 

9 = 0 a sin — t (13) 

T 

where 8 a = j-£. By specifying V t and r, one can 
obtain the amplitude of the angular motion, 9 a , for a 
certain value of the cone radius, r c . Next, we present 
the results for different values of V s and r of the 
periodic rotatory oscillation. 

1. Rotatory Oscillation; ^=0.06, r=7.2, 8 a =45°: 

The corresponding frequency of this motion is 
0.873. The results of this case are given in Fig. 10. 
The period of the Cy response is observed to be 
7.2 which is the same as that of the motion. The 
Cy curve oscillates between -0.00077 and -0.00054 
with a mean value of -0.000655, which is between 
the mean values of CCW and CW spinning cases of 
Figs. 2 and 4. Hence, these values of Vs , T and 9 a 
for the rotatory oscillation do not reduce the mean 
value of Cy. 

2. Rotatory Oscillation; V s = 0.2, r = 4.3, 9 a = 90°: 

The corresponding frequency of this motion is 
1.461. The results of this case are given in Fig. 11. 
The period of the Cy response is observed to be 4.3 
which is the same as that of the motion. The Cy curve 
oscillates between -0.00105 and -0.0002 with a mean 
value of -0.000625. This mean value of the Cy is 
higher than that of the CCW spinning case of Fig. 6 
and is slightly higher than that of the CW spinning 
case of Fig. 8. Hence, these values of Vj.r and 6 a 
for the rotatory oscillation do not reduce the mean 
value of Cy. 

3. Rotatory Oscillation; V s = 0.5, 

7 = 7.2, 9 a = 375°: 

The corresponding frequency of this motion is 
0.873. which is the same as that of the case of Fig. 10. 
However, the amplitudes of the surface velocity and 
angular motion are one order of magnitude higher 
than those of the case of Fig. 1 0. The results of this 
case are given in Fig. 12. Although the Cy response 
is periodic with the same period as that of the motion, 
there are several peaks within each period. The Cy 
changes sign from positive to negative and the mean 
value of the Cy is zero. This shows that the values of 


V,,t and 0 a for the rotatory oscillation eliminate the 
Cy on the average. It should be emphasized here that 
both the amplitudes of surface velocity and angular 
motion are one order of magnitude higher than that 
of the case of Fig. 10. 

4. V s = 0.5, r = 43, 9 a = 225°: 

The corresponding frequency of this motion is 
1.461, which is the same as that of the case of Fig. 1 1 . 
However, the amplitudes of the surface velocity and 
angular motion are 2.5 times as those of the case of 
Fig. 11. The results of this case are given in Fig. 13. 
The Cy response is periodic with the same period 
as that of the motion, but with several peaks within 
each period. Here again the Cy changes sign from 
positive to negative and the mean value of the Cy 
is -0.00015, which is better than that of any of the 
uniform spinning cases or the rotatory oscillations of 
Figs. 10 and 11. However, it is higher than that of 
the previous case. The only differences between this 
case and the previous case is the period of oscillation 
and the amplitude of angular motion. Although the 
magnitude of the mean Cy is higher than that of the 
case of Fig. 12, the peak values of Cy of the present 
case are substantially lower than those of the previous 
case. It seems that the best Cy response (zero mean 
and small amplitude) can be achieved by using the 
higher 9 a and the lower r of the present case and 
the previous case. Optimal control laws should be 
developed to effectively investigate this problem. 

CONCLUDING REMARKS 

In the present study, the effectiveness of uniform 
spinning and rotary oscillation as active control meth- 
ods for alleviating the flow asymmetry and the side 
force has been investigated computationally. It has 
been shown that a large value of uniform CCW spin- 
ning rate is very effective in substantially reducing the 
side force on the average for the given initial case of 
asymmetric flow. The CCW spinning increases flow 
separation on the left side and delays it on the right 
side, which produces equal positive and negative side 
forces within each cycle of the side-force response. 
The rotatory oscillation with large surface-velocity 
amplitude, large angular-motion amplitude and small 
period of oscillation (high frequency) is much more 
effective than the uniform CCW spinning for the same 
surface velocity because it does not only eliminate 
the mean side force but it also reduces the amplitude 
of the side force substantially. Moreover, the effec- 
tiveness of the rotatory oscillation control does not 
require certain initial shape of the vortex asymmetry. 
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Work is currently underway to study the effectiveness 
of the rotatory oscillation control on asymmetric flow 
cases with periodic vortex shedding at high angles of 
attack. Optimal control laws should be developed for 
investigating this problem. 
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(CCW spinning), a = 20°, Moo = 1*8, Re = 10 5 i T = 2.749. 
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(CW spinning), a — 20°. Moo = 1*8. Re “ 10°. r = 2.749. 
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ABSTRACT 

The effects of freestream Mach number and angle of 
a narir on the leading-edge vortex breakdown due to the 
terminating shock on a 65— degree, sharp-edged, cropped 
delta wing are investigated computationally. The compu- 
tational investigation uses the time-accurate solution of 
the laminar, unsteady, compressible, full Navier-Stokes 
equations with the implicit, upwind, flux-difference split- 
ting, finite-volume scheme. A fine O-H grid consisting of 
125x85x84 points in the wrap-around, normal and axial 
directions, respectively, is used for all the flow cases. 
Keeping the Reynolds number fixed at 3.23 xlO 6 , the 
Mach number is varied from 0.85 to 0.9 and the angle 
of attack is varied from 20° to 24°. The results show that 
at 20° angle of attack, the increase of the Mach number 
from 0.85 to 0.9 results in moving the location of the ter- 
minating shock downstream. The results also show that 
at 0.85 Mach number, the increase of the angle of at- 
tack from 20° to 24° results in moving the location of 
the terminating shock upstream. The results are in good 
agreement with the experimental data. 

INTRODUCTION 

The literature shows that vortical flows around delta 
wings in the low-speed regime have received a substantial 
volume of experimental 1 ' 4 and computational 5-9 research 
work. In the high angle of attack range, vortical flows in 
the low-speed regime are characterized with three types 
of boundary-layer separation, namely; primary, secondary 
and tertiary separations. The primary separated flow rolls 
up into a strong primary vortex core which produces a 
strong suction-pressure peak on the wing surface. The 
spanwise adverse-pressure gradient of the primary vortex 
causes the spanwise, outboard-moving, boundary-layer 
flow to separate forming a secondary vortex with opposite 
sense of rotation to and smaller strength than that of the 
primary vortex. The spanwise adverse-pressure gradient 
of the secondary vortex causes the spanwise, inboard- 
moving, boundary-layer flow to separate forming a ter- 
tiary vortex with same sense of rotation as and substan- 
tially small strength than that of the primary vortex. The 
spanwise surface-pressure curves are characterized with 
three suction-pressure peaks which varies in strength cor- 
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responding to the locations of the primary, secondary and 
tertiary vortices. When the angle of attack reaches a crit- 
ical value, the axial-pressure gradient and the high swirl 
ratio of the primary vortex produce a stagnation point 
along the path line of the primary-vortex core, and vor- 
tex breakdown of the primary core develops. Depending 
on the swirl ratio, axial pressure gradient and Reynolds 
number, the primary-core vortex -breakdown mode might 
be a bubble type, a spiral type or a bubble-spiral type. 

As the freestream Mach number increases, the vortical 
flow around the delta wing changes substantially due to 
the compressibility effects. In the supersonic flow regime, 
shock waves appear beneath or above the primary vor- 
tex, depending on the freestream normal Mach number 
and normal angle of attack. Experimental data 10, 11 and 
the computational results 12-14 have shown these types of 
vortical-flow structures. The foot print of these shock 
waves runs along a ray line from the wing vertex. If 
the shock wave is beneath the primary vortex, it interacts 
with the spanwise, outboard-moving, boundary-layer flow 
and causes, in addition to the adverse pressure gradient 
produced by the primary vortex, secondary-flow separa- 
tion. If the shock wave is above the primary vortex, it 
flattens the primary vortex and the spanwise surface pres- 
sure curve. Comparison of the surface pressure distribu- 
tion over a delta wing in low-speed and supersonic-speed 
regimes, shows that the suction-pressure peak correspond- 
ing to the primary vortex is lower for the supersonic flow 
than that for the low-speed flow. 

In the transonic-flow regime, research work on vor- 
tical flows around delta wings was given adequate at- 
tention only recently. Understanding the steady and un- 
steady, transonic, vortical-flow structures around delta 
wings in the moderate-high angle of attack range is im- 
portant for increasing the performance quality of the new 
generation of supermaneuver aircraft (e.g. YF22). Recent 
experimental measurements of transonic flows around a 
65° cropped delta wing 15-21 show that a complex shock- 
wave system appears over the upper wing surface. The 
shock-wave system consists of a ray shock wave beneath 
the leading-edge primary vortex and a transverse, time- 
dependent 16 , normal-shock wave (known as a terminat- 
ing shock) which runs from the plane of symmetry to the 
wing leading edge. The terminating shock wave interacts 
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with the primary-vortex core causing it to breakdown at 
an angle of attack as low as 18°. Such a critical angle 
of attack is substantially smaller than the critical angle 
of attack of vortex breakdown in the low-speed regime. 
Reference 21 contains extensive flow measurements for 
the 65° cropped delta wing with and without leading- 
edge extension (LEX). A complete reconstruction of the 
three-dimensional flow field at and behind the terminating 
shock was not possible experimentally. 

Computational simulations for transonic delta-wing 
flows have been developed on a very limited scale by 
using the Euler equations 20 * 22 and the thin-layer Navier- 
Stokes equations 23 . The Euler-equations solutions were 
not capable of fully resolving the flow in the terminating 
shock region and the thin-layer Navier-Stokes-equations 
solutions did not address that region. In Ref. 24 by the 
present authors, the laminar, unsteady, compressible, full 
Navier-Stokes equations are integrated time accurately us- 
ing the implicit, upwind, flux-difference splitting, finite- 
volume scheme to study and construct the flow field 
structure of a transonic flow around a 65° sharp-edged, 
cropped-delta wing at 20° angle to attack, 0.85 Mach 
number and 3.23 xlO 6 Reynolds number. A fine O-H 
grid consisting of 125x85x84 points in the wrap-around, 
normal and axial directions, respectively, is used for the 
computational solution. A A-shock system, which con- 
sists of a ray shock under the primary vortex core and 
a transverse terminating shock, has been captured. Be- 
hind the terminating shock, the leading-edge vortex core 
breaks down into a two-bubble cell type. The terminat- 
ing shock and the vortex breakdown region behind it are 
time dependent and appear to be oscillatory. The flow 
field ahead of the terminating shock is steady and in- 
cludes a supersonic pocket which is surrounded by the 
ray shock and the terminating shock. The flow inside 
the pocket does not change due to changes in the flow 
downstream. This is consistent with the fact that the su- 
personic pocket along with the terminating shock do not 
allow disturbances to propagate upstream. These results 
have been validated using the available experimental data 
and they are in good agreement. This work gives a com- 
plete construction of the flow field over the wing surface 
and in particular the structure of the flow at the terminat- 
ing shock and behind it. 

In this paper, a parametric study is carried out to in- 
vestigate the effects of frees tream Mach number and an- 
gle of attack on the terminating shock and the leading- 
edge, primary-vortex breakdown for the same 65° sharp- 
edged, cropped delta wing. The computational investiga- 
tion uses the same equations, computational scheme and 
grid of Ref. 24. Keeping the Reynolds number fixed at 
3.23 x 10 6 , the Mach number is changed from 0.85 to 0.9 
while the angle of attack is fixed at 20°, and the angle 
of attack is changed from 20° to 24° while the Mach 
number is fixed at 0.85. 


HIGHLIGHTS OF FORMULATION 
AND COMPUTATIONAL SCHEME 

The conservative form of the dimensionless, unsteady, 
compressible, full Navier-Stokes equations is used for the 
formulation of the problem. The equations are written in 
terms of the time-independent, body-conformed coordi- 
nates and (Ref. 25). 

The implicit, upwind, flux -difference splitting, finite- 
volume scheme is used to solve the unsteady, compress- 
ible, full Navier-Stokes equations. The scheme uses the 
flux-difference splitting scheme of Roe which is based 
on the solution of the approximate one-dimensional, Rie- 
mann problem. In the Roe scheme, the inviscid flux dif- 
ference at the interface of computational cells is split into 
two parts; left and right flux differences. The splitting is 
accomplished according to the signs of the eigenvalues of 
the Roe averaged-Jacobian matrix of the inviscid fluxes 
at the cell interface. The smooth flux limiter is used to 
eliminate oscillations at locations of large flow gradients. 
The viscous- and heat-flux terms are linearized in time 
and the cross-derivative terms are neglected in the im- 
plicit operator and retained in the explicit terms. The vis- 
cous terms are differenced using a second-order accurate 
central differencing. The resulting difference equation is 
approximately factored and is solved in three sweeps in 
the and directions. The computational scheme 
is coded in the computer program “FTNS3D” which is a 
modified version of the CFL3D-code. 

COMPUTATIONAL RESULTS 

A 65° swept-back, sharp-edged, cropped delta wing 
of zero thickness is considered for the computational so- 
lutions. The cropping ratio (tip length/root-chord length) 
is 0.15. An O-H grid of 125 x 85 x 84 in the wrap-around, 
normal and axial directions, respectively, is used. The 
computational domain extends two-chord length forward 
and five-chord length backward from the wing trailing 
edge. The radius of the computational domain is four- 
chord length. The minimum grid size normal to the wing 
surface is 5x10^ from the leading edge to the plane of 
symmetry. Figure 1 shows a three-dimensional shape of 
the grid and a cross-flow plane. 

Time-accurate integration of the laminar, unsteady, 
compressible, full Navier-Stokes equations has been car- 
ried out with At = 0.0002. Three flow conditions are 
used to study the effect of increasing the Mach number 
while the angle of attack is kept constant and the effect 
of increasing the angle of attack while the Mach num- 
ber is kept constant In all the three cases, the Reynolds 
number, R*, is 3.23 x 10 6 based on the root-chord length. 

Case I (Moo = 0.85, a = 20°) 

For this case, the freestream Mach number, Af*,* and 
angle of attack, a, are 0.85 and 20°, respectively. Figure 2 
shows a comparison of the computed, span wise, surface- 
pressure coefficient (C p ) at different chord stations (x = 
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0.3, 0.6 and 0.8) with the experimental data of Erickson 21 
(Re = 3.23 xlO 6 ) and Hartmann 17 (Re = 2.38 xlO 6 and 
4.57 xlO 6 ). The computational results show the correct 
location and level of the suction-pressure peak corre- 
sponding to the primary vortex in comparison with the 
experimental data. They also show a smaller suction- 
pressure peak corresponding to the secondary vortex. The 
computational results are in fair to good agreement with 
the experimental data. For the chord station x = 0.9, the 
Cp -curve shows a rapid increase in the pressure coeffi- 
cient (a decrease in the suction-pressure coefficient). For 
example, the suction-pressure-peak coefficient increases 
from a value of -1.4 at x = 0.8 to a value of -1.15 at x = 
0.9. Figure 3 shows the total-Mach contours and stream- 
lines at the chord stations of x = 0.60, 0.90 and 0.97. 
At x = 0.60, the Mach contours show an oblique shock 
beneath the primary vortex and a subsonic, separated re- 
gion to its right. The streamlines show a secondary sepa- 
rated flow and the corresponding secondary vortex. This 
separation is due to the shock interaction with the sur- 
face boundary-layer flow and is also due to the adverse, 
spanwise pressure gradient created by the primary vor- 
tex. At x = 0.90, the shock beneath the primary vortex 
becomes weak and the primary -vortex size increases. At 
x = 0.97, the shock beneath the primary vortex disap- 
pears and the primary vortex diffuses and reduces to a 
repelling focus, as shown by the streamlines. The details 
of the flow structure at x = 0.90 and 0.97 in addition 
to the spanwise, pressure-distribution curve at x = 0.90 
clearly indicate that the primary vortex is experiencing a 
vortex breakdown due to a transverse shock (terminating 
shock) which is located between x = 0.80 and x = 0.90. 

Figure 4 shows the static pressure contours on the 
wing and symmetry planes. The contours clearly show 
the location, shape and strength of the terminating shock. 
A substantial supersonic pocket which is bounded by the 
terminating shock and the ray shocks (shocks beneath 
the primary-vortex cores) is observed on the wing plane. 
The terminating shock is located at x = 0.83 at the 
plane of symmetry, which is in good agreement with the 
experimental data 21 , where the shock is located at x = 0.84 
at the plane of symmetry. Figure 5 shows the position of 
ray lines from the wing vertex (which are marked by the 
letters A-H) and the static-pressure variation along these 
lines. The static-pressure curves give several points to 
generate the foot-print line of the terminating shock. The 
terminating shock is found to extend from the plane of 
symmetry to the wing leading edge. It reaches its highest 
strength at the location of the primary vortex (lines E-G). 
Figure 6 shows the total-Mach contours and streamlines 
on a vertical ray plane at the 0.68 spanwise location which 
passes through the vortex breakdown. Blow-ups of the 
velocity vectors and streamlines on this ray plane are also 
shown in Fig. 6. The streamlines conclusively show a 
two-bubble cell vortex breakdown. It is a typical three- 
dimensional vortex breakdown mode which consists of 
an attracting saddle point (front), a repelling saddle point 
(rear), an attracting focus (top), and a repelling focus 


(bottom). Such a breakdown mode is similar to the one 
which was captured for an isolated supersonic vortex in 
an unbounded domain in Refs. 26 and 27. The location 
of the attracting saddle point is at 0.97 along the ray line 
which corresponds to a location of 0.87 along the axial 
direction. The Mach contours show that the front surface 
of the vortex-breakdown bubbles is enclosed by a hemi- 
spherical shape-like shock surface. In Fig. 18, the details 
of the flow structure on the wing and symmetry planes 
are shown. 

Having established the flow structure of this case, the 
Mach number is increased to 0.9 while the angle of attack 
is kept fixed at 20°. 

Case II (Moo = 0.90, a = 20°) 

The results of this case are given in Figs. 7-11 and 
19. Figure 7 shows the computational spanwise, surface- 
pressure coefficient at different chord stations along with 
the experimental data of Erickson 21 . The computational 
results are in good agreement with the experimental data 
at x = 0.3 and 0.6. At x = 0.8, the computational re- 
sults underestimates the pressure coefficient of the exper- 
imental data. The locations of the primary and secondary 
vortex cores are in good agreement with those of the ex- 
perimental data. It is noticed that the levels of C p for 
the present case are lower than those of Case I (Fig. 2). 
Again, the pressure level decreases rapidly at x = 0.90. 
Figure 8 shows the total-Mach contours and streamlines 
in cross-flow planes at x = 0.60, 0.90 and 0.97. The shock 
beneath the primary vortex is observed in the Figures at 
x = 0.60 and x = 0.90. For x = 0.90, the shock beneath 
the primary vortex is still strong in comparison with that 
of Case I (Fig. 3). At x = 0.97, the repelling focus is 
observed indicating that vortex breakdown has occurred. 
Figure 9 shows that the terminating shock in the cross- 
flow plane is located at x = 0.93 within the boundary - 
layer, which is in good comparison with the experimen- 
tally measured shock of Ref. 21, where it is located at x 
= 0.95. The static -pressure contours on the wing plane 
show that the terminating shock for Case II (Fig. 9) is 
closer to the trailing edge that of Case I (Fig. 4). It should 
be noted here that the terminating -shock location in the 
outer flow is ahead of its location in the boundary-layer 
flow. The static -pressure variations along the ray lines of 
Fig. 10 clearly show that the terminating-shock foot print 
is located between x = 0.925 and x = 0.95, and that it 
extends from the plane of symmetry to the wing leading 
edge. Figure 11 shows the Mach contours and stream- 
lines on a vertical ray plane passing through the vortex 
breakdown. It is noticed that the vortex breakdown shape 
is different from and smaller than that of Case I (Fig. 6). 
The attracting saddle point, attracting focus and repelling 
saddle point are clearly observed. The repelling focus 
is very small. This indicates that the terminating shock 
becomes smaller in strength than that of Case I. Figure 
19 shows the details of this flow case on the wing and 
symmetry planes. 
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It is concluded that as the freestream Mach number 
increases slightly from 0.85 to 0.9, the terminating shock 
strength decreases and its location moves downstream 
from x = 0.84 to x = 0.93. Moreover, the surface pressure 
levels become smaller than those of Case I. 

Next, the Mach number is kept fixed at 0.85 and the 
angle of attack is increased to 24°. 

Case in (Moo = 0.85, a = 24°) 

The results of this case are given in Figs. 12-17 and 
20. The computational surface-pressure result at x = 0.3 
(Fig. 12) is in good agreement with the experimental data 
of Erickson 21 . However, the computational results, at x 
= 0.6 and 0.8 are either overpredicting or underpredicting 
the experimental data. Figures 13, 14 and 15 show that 
the terminating shock moves upstream to x = 0.753 in the 
boundary-layer flow at the plane of symmetry. This is in 
good agreement with the experimental data of Ref. 21, 
where the shock is located at x = 0.75 in the boundary 
layer flow. The terminating-shock location in the outer 
flow is ahead of its location in the boundary layer. Figure 
16 shows that the vortex-breakdown region is larger than 
those of Cases I and II. Moreover, the attracting and 
repelling foci are smaller than those of Case I. Figure 20 
shows die details of this case on the wing and symmetry 
planes. 

Thus, it is seen that as the angle of attack increases 
from 20° to 24° while the Mach number is kept fixed 
at 0.85, the terminating shock moves upstream and the 
vortex-breakdown region becomes large. Moreover, the 
surface pressure levels become larger than those of Case I. 

The computational results show that the flow at the 
terminating shock and behind it is time dependent and it 
indicates oscillatory motion (The computations have not 
been carried out beyond t = 6.0 or 30,000 time steps 
with At = 0.0002). In Fig. 17, we show snapshots of the 
streamlines and their blow-ups on a ray plane passing 
through the vortex-breakdown region. The snapshots 
are shown at t = 4.22, 5.16 and 5.52. It is clearly 
seen that the vortex breakdown moves upstream showing 
different modes. In the same time, the terminating shock 
is also moving upstream and slows down to reverse its 
direction of motion. This is in complete agreement with 
the experimental observations of Bannik and Houtmann 16 . 

CONCLUDING REMARKS 

The laminar, unsteady, compressible, full Navier- 
Stokes equations are integrated time accurately using the 
implicit, upwind, flux-difference splitting finite-volume 
scheme to study the transonic flow field around a 65° 
sharp-edged, cropped delta wing. First, the flow field has 
been constructed for a Reynolds number of 3.23 xlO 6 , 
a Mach number of 0.85 and an angle of attack of 20° 
(Case I). A A-shock system consisting of a ray shock be- 
neath the primary vortex core and a transverse terminating 
shock has been captured. Behind the terminating shock, 


the leading-edge vortex core breaks down. Keeping the 
Reynolds number constant and the angle of attack fixed 
at 20° , the Mach number is increased to 0.90. The results 
of this case (Case II) show that the terminating shock 
moves downstream and the vortex -breakdown region be- 
comes smaller than that of Case I. Keeping the Reynolds 
number constant and the Mach number fixed at 0.85, the 
angle of attack is increased to 20°. The results of this 
case (Case III) show that the terminating shock moves up- 
stream and the vortex-breakdown region becomes larger 
than that of Case I. The computational results are in good 
agreement with the experimental data. However, it must 
be emphasized that the flow at the terminating shock and 
behind it is time dependent while the flow ahead of the 
terminating shock is steady. The present paper shows the 
structure of the flow field behind the terminating shock 
for the first time. 
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Fig. 7 Comparison of the computed and experimental spanwise, surface-pressure 
coefficient at different chord stations; Moo = 0.90, a = 20°. 
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Fig. 12 Comparison of the computed and experimental spanwise, surface-pressure 
coefficient at different chord stations; Moo = 0.85, tv = 24°. 
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Fig. 18 Surface-pressure and Mach contours and particle trace on wing and 
symmetry planes; .l/ x = 0.85, n = 20°. 
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Fig. 19 Surface-pressure and Mach contours and particle trace on wing and 
symmetry planes: _U X = 0.90, n = 20° 




Supersonic Vortex Breakdown on a Delta Wing 
M = 0.85, Re = 3,230,000 and AQA = 24* 


















Surface Pressure 


' ..Ai" 


0.11 0.52 0.93 

M 0.85, Re - 3,230,000 and AO A ^24 


*K 


a**®5rs»v; 

: • • ■; . ;; 


'3 

ii: .•:■• : * ?. .. ; ®“® 

a : *> 

i r* . <; 


, >; \. { 


?: '■ i- 




Mach Contours 

ii 

0.75 1.02 130 


'■ git!, 

. w- 




’i# v n. 


Fig. 20 Surface-pressure and Mach contours and particle trace on wing and 
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ABSTRACT 

Transonic flow over a 65-degree swept-back, sharp- 
edged, cropped delta wing is investigated computationally 
nsing the time-accurate solution of the unsteady, com- 
pressible, full Navier-Stokes equations with an implicit, 
upwind, flux-difference splitting, finite-volume scheme. 
Coarse and fine O-H grids are used to obtain the so- 
lution. The grid consists of 125 x 85 x84 points in the 
wrap-around, normal and axial directions, respectively. 
The results are presented for an angle of attack of 20°, 
March number of 0.85 and Reynolds number of 3.23 x 10 6 
(based on the wing chord length). With the fine grid, the 
results show that a system of shocks has been captured 
over the upper wing surface, and that the leading-edge 
vortex core experiences an unsteady, supersonic vortex 
breakdown after passing through a spanwise shock (ter- 
minating shock) near the wing trailing edge. The com- 
puted results at a certain time are in good agreement with 
the experimental data. Topological aspects of the vortex- 
breakdown flowfield are also presented and discussed. 

INTRODUCTION 

At sufficiently high angles of attack, vortex break- 
down for incompressible flows around delta wings has 
been observed along the leading-edge primary vortex 
cores. Two distinct forms of vortex breakdown have been 
documented experimentally 1 . The first form is the bubble 
type and the second form is the spiral type. The bub- 
ble type shows an almost axisymmetric sudden swelling 
of the vortex core into a bubble, while the spiral type 
shows an asymmetric, spiral, vortex filament followed by 
a rapidly spreading turbulent flow. Both types are charac- 
terized by an axial stagnation point and a limited region 
of reversed axial flow. Much of our knowledge of in- 
compressible vortex breakdown has been obtained from 
experimental studies of pipe flows where both types of 
breakdown and other types as well were generated and 
documented 2-4 . 

The major effort of computational study of vor- 
tex breakdown flows has also— been -focused on iso- 
lated swirling flows. For incompressible flows, quasi- 
axisymmetric, bubble-type, vortex-breakdown flows were 
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fpmpnturf using the Navier-Stokes equations 5- *. Three- 
dimensional bubble and spiral vortex-breakdown flows 
were t»i»» computed for isolated swirling flows using 
the three-dimensional Navier-Stokes equations in the 
vorticity-velocity form or the primitive variables form’ -11 . 

In ter action between a longitudinal vortex and a trans- 
verse shock wave occurs in several flow applications 
which include transonic and supersonic flows over a delta 
wing tv s stroke- wing configuration at moderate to high 
an gle s of a supersonic inlet ingesting a vortex, 

and a supersonic combustor where fuel is injected in a 
swirling jet to enhance fuel-air mixing 15-14 . For delta 
wings and stroke- wing configurations, vortex breakdown 
is an undesirable phenomenon since it produces wing 
staff Therefore, its occurrence needs to be delayed with 
passive or active control methods in order to increase 
the wing performance at large angles of attack. For 
a supersonic combustor, vortex breakdown is desirable 
since it enhances mixing of air and fuel and stabilizes the 
flame 15,16 . Therefore its occurrence needs to be enhanced 
and controlled. 

For supersonic flows, quasi-ax is ymmetric bubble-type 
vortex-breakdown 17-19 and three-dimensional bubble-type 
and spiral-type vortex breakdown 20 for isolated swirling 
flows have been recently computed by the present au- 
thors. Using compatible, inlet boundary conditions, the 
time-accurate solutions of the unsteady, compressible, full 
Navier-Stokes equations were obtained to study the ef- 
fects of Reynolds number, Mach number, swirl ratio, 
type of exit-boundary conditions and grid fineness and 
distribution on the vortex-breakdown modes for internal 
and external flows. Several modes of vortex breakdown 
which include transient single-bubble, transient multi- 
bubble, periodic multi-frequency multi-bubble, quasi- 
steady two-bubble cell and spiral-type vortex break- 
downs have been obtained 21 . For three-dimensional 
vortex-breakdown flows in a swirling, supersonic jet 
flow, topological aspects of the critical points in the 
vortex-breakdown region have been studied and com- 
pared with the available experimental incompressible 
vortex-breakdown topology. 

Recent experimental measurements 22-26 of transonic 
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flows around a 65° swept-back, cropped delta wing show 
that shock wave formation is likely to occur underneath 
the leading-edge primary vortex core. In cross-flow 
planes perpendicular to the wing, the cross-flow beneat h 
the primary vortex reaches supersonic speeds and a cross- 
flow shock develops beneath the primary vortex similar 
to the supersonic flow in a convergent-divergent nozzle. 
These measurements also show that a transverse shock 
‘Terminating shock” which might cause primary-vortex- 
core breakdown could develop in an analogous manner 
to the shock that terminates the two-dimensional super- 
sonic pocket on an airfoil. A complete reconstruction 
of the three-dimensional flow field on the delta wing in 
this region was not possible experimentally 33 ' 36 . Com- 
putational simulations for transonic delta-wing flows have 
been developed by using the Euler equations 37, * and 
the thin-layer, Navier-Stokes equations 29 . The Euler- 
equations solutions were not capable of fully resolving 
the flow in the terminating shock region and the thin- 
layer, Navier-Stokes solution did not address that region. 

In the present paper, we consider the transonic flow 
around a 65° sharp-edged, cropped delta wing at an angle 
of attack of 20°, a Mach number of 0.85 and a Reynolds 
number of 3.23 x 10 6 . The purpose of the present numer- 
ical simulation and study is to construct the flow field 
o ver the wing with particular emphasis of the vortex 
core-terminating shock interaction region. The laminar, 
unsteady, compressible, full Navier-Stokes equations are 
solved accurately in time with an implicit, flux-difference 
splitting, finite-volume scheme. The computations are 
carried out with time-accurate stepping on two O-H grids; 
a coarse grid and a fine grid. Both grids consist of 
125x85x84 points in the wrap-around, normal and axial 
directions, respectively. The main difference between the 
coarse and fine grids is the distribution of the grid points 
normal to the wing surface within the thin viscous layer 
(to be discussed later on). 

HIGHLIGHTS OF FORMULATION 

AND COMPUTATIONAL SCHEME 

The conservative form of the dimensionless, unsteady, 
compressible, full Navier-Stokes equations is used for the 
formulation of the problem. The equations are written in 
terms of the time-independent, body-conformed coordi- 
nates ,£ 2 and (Ref. 18). 

The implicit, upwind, flux-difference splitting, finite- 
volume scheme is used to solve the unsteady, compress- 
ible, full Navier-Stokes equations. The scheme uses the 
flux-difference splitting scheme of Roe which is based 
on the solution of the approximate one-dimensional, Rie- 
mann problem. In the Roe scheme, the inviscid flux dif- 
ference at the interface of computational cells is split into 
two pans; left and right flux differences. The splitting is 
accomplished according to the signs of the eigenvalues of 
the Roe averaged-Jacobian matrix of the inviscid fluxes 
at the cell interface. The smooth flux limiter is used to 
eliminate oscillations at locations of large flow gradients. 


The viscous- and beat-flux terms are linearized in time 
and the cross-derivative terms are neglected in the im- 
plicit operator and retained in the explicit terms. The vis- 
cous terms are differenced using a second-order accurate 
central differencing. The resulting difference equation is 
approximately factored and is solved in three sweeps in 
the (M 3 and { 3 directions. The computational scheme 
is coded in die computer program “FTNS3D” which is a 
modified version of the CFL3D-code. 

COMPUTATIONAL RESULTS 

A 65° swept-back, sharp-edged, cropped delta wing 
with zero thirknftsj; is considered for the computational 
solutions. The cropping ratio (tip length/root-chord 
length) is 0.15. The wing angle of attack is 20°, and the 
freestream Mach number and Reynolds number (based 
on die root-chord length) are 0.85 and 3.23 xlO 6 , re- 
spectively. The reason behind the present, selected flow 
cnnrii riniw is because of the uncertainty of the existing 
experimental data 23-26 about the structure of the down- 
stream flow field of the leading-edge vortex core. The 
experimental data shows that a supersonic flow region 
appears on the upper wing surface near the plane of sym- 
metry. This flow region is terminated by a transverse 
shock (known as a terminating shock) in a similar way to 
the shock that terminates a supersonic pocket on a super- 
critical airfoil 23 . 

Grid: 

An O-H grid of 125x85x84 in the wraparound, 
normal and axial directions, respectively, is used for the 
computational simulation. The computational domain 
extends two-chord length forward and five-chord length 
backward from the wing trailing edge. The radius of the 
computational domain is four-chord length. Two grids 
have been constructed using the same number of grid 
points. The first is called the coarse grid and the second 
is called the fine grid. For the coarse grid, the grid points 
in the cross flow planes have been distributed using a 
Joukowski transformation which produces a minimum 
grid size, normal to the wing surface, that varies from 
5x10^ at the leading edge to 3xl(T 2 at the plane of 
symmetry. For the fine grid, the elliptical grid lines in 
the cross-flow planes have been constructed such that the 
minimum grid size normal to the wing surface, stays 
constant at SxlO -4 from the leading edge to the plane 
of symmetry. Figures 1 and 2 show three-dimensional 
shape of the coarse and fine grids and a cross-flow plane 
along with its blow-ups. 

Time-accurate integration of the laminar, unsteady, 
compressible, full Navier-Stokes equations has been car- 
ried out with At = 0.001 for the coarse grid and At * 
0.0002 for the fine grid. The results showed that the 
leading-edge vortex core passes through a terminating 
shock which causes the vortex core to breakdown. More- 
over, it is shown that the flow becomes unsteady behind 
the terminating shock. 
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Validation of Surface Pressure: 

Figure 3 shows a comparison of the computed, span- 
wise surface-pressure coefficient (Cp) at different chord 
stations for the fine and coarse grids with the experimen- 
tal data of Erickson 30 (R* = 3.23 xlO 6 ) and Hartmann 34 
(Re = 2.38 xlO 6 and 4.57 xlO 6 ). The computed results 
are selected at t = 3.6. Obviously, the coarse-grid Cp- 
curves do not show the suction-pressure peak correspond- 
ing to the secondary vortex and the correct location of the 
suction-pressure peak corresponding to the primary vor- 
tex. The coarse-grid Cp -curves are similar to those of the 
Euler-equations solution. Therefore, they are discarded 
in this paper. The fine-grid C p -curves show the correct 
location of the suction-pressure peak corresponding to the 
primary vortex and the suction-pressure peak correspond- 
ing to the secondary vortex. The fine-grid Cp-curves at 
x = 0.3, 0.6 and 0.8 arc in fair to good agreement with 
the experimental data. For x = 0.9, the fine-grid Cp- 
curve shows a substantial, rapid increase in the pressure 
coefficient (a decrease in the suction pressure). Figure 
4 shows the total-Mach contours and the streamlines in 
cross-flow planes at the chord stations of x = 0.3, 0.6, 0.8, 
0.9, 0.97 and 1.0. At x = 0.3, 0.6 and 0.8, the total-Mach 
contours show an oblique shock under the primary vortex 
and a small subsonic region to the right of the shock. The 
streamlines show the secondary separation to the right of 
the shock. This separation is due to the shock interac- 
tion with the surface boundary-layer flow and is also due 
to the adverse, spanwise pressure gradient created by the 
primary vortex. At x = 0.9, the shock under the primary 
vortex becomes weak as observed in the total-Mach con- 
tours and the primary-vortex size increases. At x = 0.97, 
the shock under the primary vortex disappears and the 
primary vortex diffuses and reduces to a repelling focus 
as shown by the streamlines. At x = 1.0, the repelling 
focus becomes a repelling line. The details of the flow 
structure shown at x = 0.9, 0.97 and 1.0 indicate that 
the primary vortex is going through a breakdown mode 
which is caused by a transverse shock (terminating shock) 
between x = 0.8 and x = 0.9. 

Terminating Shock: 

To show that a terminating, transverse shock exists 
and has been captured computationally, the static -pressure 
contours and total-Mach-contours on two planes are com- 
puted and displayed in Fig. 5. In this figure, the static- 
pressure contours are shown on the wing surface and 
the plane of symmetry, and the total-Mach contours are 
shown on the third plane (k = 3) above the wing (in the 
viscous layer) and on the plane of symmetry. The plane 
of symmetry contours clearly show the location, shape 
and strength of the terminating shock. Moreover, the 
Mach contours show that a substantial supersonic pocket 
(bounded by the sonic line and terminating shock) ex- 
tends from the wing vertex to the shock location of x = 
0.83, which is in good agreement with the experimental 
data 30 , where the shock is located at x = 0.84. The com- 


puted results show that the shock is a normal shock with 
a height of 0.4 which is equal to one-half the wing span. 
In the spanwise direction, the shock foot print (shown on 
the Mach contours at k « 3) extends beyond the primary- 
vortex location. A A-type shape of the shock-system foot 
print, which on one side of die wing, consists of the ter- 
minating shock and the shock under the primary vortex 
that runs along a ray plane from the wing vertex, is seen 
on the Mach contours at k = 3. 

Figure 6 shows the position of the ray lines from the 
wing vertex (which arc marked by the letters A-H) and 
the static-pressure curves along these lines. The static- 
pressure curves show the spanwise locations of several 
points on the foot-print line of the terminating shock. The 
terminating shock is clearly seen to run in the spanwise 
direction from the plane of symmetry to the wing leading 
edge. It reaches its highest strength from the location of 
the primary vortex to the wing leading edge (from line 
E to line H). 

Vortex-Breakdown Structure: 

Having established the shock system that consists of 
the shock under the primary vortex and the terminating 
shock, the focus is directed on the structure of the flow be- 
hind the terminating shock. In Fig. 7, we show the total- 
Mach contours and streamlines on a ray plane at the 0.658 
spanwise location, which passes through the leading-edge 
vortex core. Blow-ups of the velocity vectors and stream- 
lines on this vertical plane are also shown in Fig. 7. The 
streamlines figures clearly show a two-bubble cell vortex 
breakdown. This is a typical three-dimensional vortex- 
breakdown mode which consists of an attracting saddle 
point (front) a repelling saddle point (rear), an attracting 
focus (top) and a repelling focus (bottom). Such a break- 
down mode is similar to the one which was captured for 
an isolate^ supersonic vortex in an unbounded domain 
in Refs. 20 and 21. The location of the attracting sad- 
dle point is at 0.97 along the ray line, which corresponds 
to 0.87 along the axial direction. The attracting focus 
point is characterized with spiralling-in streamlines and 
the repelling focus point is characterized with spiralling- 
out streamlines. The Mach contours show that the front 
surface of the vortex-breakdown bubbles is enclosed by a 
hemi-spherical shape-like shock surface. Figures 12 and 
13 show details of the flow structure on the wing plan 
view, on the plane of symmetry and on the ray plane at 
the 0.658 spanwise location (marked as J = 16 on Fig. 13). 
These figures and discussion give a complete construction 
of the flow structure including the shock system and its 
interaction with the leading-edge vortex core which pro- 
duces vortex-breakdown of the two-bubble-cell mode. 

Unsteadiness of the Vortex-Breakdown: 

The computations have been carried out with time- 
accurate stepping beyond t = 3.6. Figures 8—11 show the 
results at t = 5,52. These results show that the terminating 
shock moves in the upstream direction and so is the 
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two bubble-cell vortex breakdown behind the terminating 
shock. Figure 8 shows that the repelling focus is at x « 

0.88 instead of x * 0.97 (Fig. 4). Figure 9 shows that 
the terminating shock in the plane of symmetry is at x ■ 
0.685 instead ofx = 0.83 (Fig. 5). The shock decreases in 
height and its thickness increases. Figure 10 shows that 
the size of the two-bubble cell vortex-breakdown region 
increases in comparison with the size at t = 3.6 (Fig. 7). 
Upstream of the terminating shock die flow stayed steady 
without any change. 

Beyond the time t * 5.52, the upstream shock motion 
stopped and the motion reversed its direction to the down- 
stream. The computations were not carried out beyond 
this instant due to its impeding cost. The unsteadiness of 
the terminating shock and the vortex-breakdown region 
behind it have also been observed experimentally by Ban- 
nik and Houtmann 23 . They also observed that the flow 
upstream of the terminating shock stayed steady without 
any change. These experimental observations undoubt- 
edly support and validate our computational results. 

CONCLUDING REMARKS 

The laminar, unsteady, compressible, full Navier- 
Stokes equations are integrated time accurately using the 
implicit, upwind, flux-difference splitting, finite-volume 
scheme to study and construct the flow field structure of 
a transonic flow around a 65° sharp-edged, cropped-delta 
wing. A A-shock system, which consists of a ray shock 
under the primary vortex core and a transverse terminat- 
ing shock, has been captured. Behind the terminating 
shock, the leading-edge vortex core breaks down into a 
two-bubble cell type. The terminating shock and the vor- 
tex breakdown region behind it is time dependent and 
appears to be oscillatory. The flow field ahead of the ter- 
minating shock stays steady without any change. This is 
consistent with the fact that the supersonic pocket along 
with the terminating shock do not allow disturbances to 
propagate upstream. The present results have been vali- 
dated using the available experimental data and they are 
in good agreement The present paper gives a complete 
construction of the flow field over the wing surface and 
in particular the structure of the flow at the terminating 
shock and behind it. 

ACKNOWLEDGEMENT 

For the first two authors, this work is supported by the 
NASA-Langley Research Center under grant No. NAG- 
1-994 along with a partial support from the AFOSR. The 
Computational Resources providetTby The“NAS Center 
at Ames and the NASA Langley Research Center are 
acknowledged and appreciated. 

REFERENCES 

1. Lamboume, N. C. and Bryer, D. W„ “Bursting of 

Leading-Edge Vortices: Some Observations and Dis- 


cussion of the Phenomenon," Aeronautical Research 
Council. RAM 3282, 1961. 

2. Sarpkaya, T., “Effect of The Adverse Pressure Gradi- 
ent on Vortex Breakdown," AIAA Journal, Vol. 12, 
No. 12, Dec. 1974, pp. 602-607. 

3. Escudier, M. P. and Zender, N., “Vortex Flow 
Regimes,” Journal of Fluid Mechanics, Vol. 115, 
1982. pp. 105-122. 

4. Leibovich, S., “Vortex Stability and Breakdown: Sur- 
vey and Extension," AIAA Journal, Vol. 22, No. 9, 
Sept 1984, pp. 1192-1206. 

5. Grabowsld, W. J. and Berger, S. A., “Solutions of 
the Navier-Stokes Equations for Vortex Breakdown,” 
Journal of Fluid Mechanics, Vol. 75, Part 3, 1976, 
pp. 525-544. 

6. Hafez, M„ Kuruvila, G. and Salas, M.D., “Numeri- 
cal Study of Vortex Breakdown,” Journal of Applied 
Numerical Mathematics, No. 2, 1987, pp. 292-302. 

7. Salas, M. D. and Kuruvila, G., “Vortex Breakdown 
Simulation: A Circumspect Study of The Steady, 
Laminar, Axi symmetric Model,” Computers and Flu- 
ids, Vol. 17. No. 1, 1989, pp. 247-262. 

8. Wu, J. C. and Hwang, S., “Computational Study of 
Vortex Breakdown in Circular Tube," AIAA 91-1820, 
June 1991. 

9. Menne, S. and Liul, C. HI., “Numerical Simulation 
of a Three-Dimensional Vortex Breakdown,” Z. Flug- 
wiss. Weltraumforsch 14, 1990, pp. 301-308. 

10. Spall, R. E., Gatski, T. B. and Ash, R. L., “The Struc- 
ture and Dynamics of Bubble-Type Vortex Break- 
down," Proc. R. Soc., London, A429, 1990, pp. 
613-637. 

11. Breuer, M. and Hanel, D., “Solution of The 3-D In- 
compressible Navier-Stokes Equations for the Simu- 
lation of Vortex Breakdown,” Eight GAMM Confer- 
ence, Delft, Netherlands, September 27-29, 1989. 

12. Krause, E., “Vortex Breakdown: Physical Issues and 
Computational Simulation,” Third International Con- 
gress of Fluid Mechanics, Cairo, Egypt, January 1990, 
Vol. 1, pp. 335-344. 

13. Delery, J., Horowitz, E., Leuchter, O. and Solignac, 
J. L., “Fundamental Studies of Vortex Flows," La 
Recherche Aerospatiale, No. 1984-2, 1984, pp. 
1-24. 

14. Metwally, O., Settles, G. and Horstman, C., “An Ex- 
perimental Study of Shock Wave/Vortex Interaction,” 
AIAA 89-0082. January 1989. 

15. Cutler, A. D. and Levey, B. S., "Vortex Breakdown 
in a Supersonic Jet,” AIAA 91-1815, June 1991. 


4 


16. Rhode, D. L., Lilley, D. B. and McLaughlin, D. K., 
"On The Prediction of Swirling Flowfields Found in 
Axisynunetric Combustor Geometries,’’ Transactions 
of ASME, Vol. 104, September 1982, pp. 378-384. 

17. KandiL O. A., Kandil, H. A. and Liu, C. H., “Compu- 
tation of Steady and Unsteady Compressible Quasi- 
Axisymmetric Vortex Row and Breakdown,” A1AA 
91-0752, January 1991. 

18. Kandil, O. A., Kandil, H. A. and Liu, C. H., “Su- 
personic Quasi- Axisymmetric Vortex Breakdown,” 
A1AA 91-331 1-CP, September 1991, pp. 851-863. 

19. Kandil, O. A., Kandil, H. A. and Liu, C. H., “Crit- 
ical Effects of Downstream Boundary Conditions on 
Vortex Breakdown,” AIAA 92-2601 CP, June 22-24, 
1992, pp. 12-25. 

20. Kandil, O. A., Kandil, H. A. and Liu, C. H., “Three- 
Dimensional Supersonic Vortex Breakdown,” AIAA 
93-0526, January 11-14, 1993. 

21. Kandil, H. A., “Navier-Stokes Simulation of Quasi- 
Axisymmetric and Three-Dimenisonal Supersonic 
Vortex-Breakdown,” Ph.D. Dissertation, Dept, of Me- 
chanical Engineering and Mechanics, Old Dominion 
University, Norfolk, VA, May 1993. 

22. Boersen, S. J. and Elsenaar, A., ‘Tests on the AFWAL 
65° Delta Wing at NLR: A Study of Vortex Row De- 
velopment Between Mach = 0.4 and 4,” Proceedings 
of Symposium on International Vortex Row Experi- 
ment on Euler Code Validation, Stockholm, Sweden, 
October 1-3, 1986, pp. 23-36. 

23. Bannik, W. J. and Houtman, E. M., “Experiments 
on the Transonic Row Over a Delta Wing at High 
Angles of Attack,” Proceedings of Symposium on 
International Vortex Row Experiment on Euler Code 


Validation, Stockholm, Sweden, October 1-3, 1986, 
pp. 37-46. 

24. Hartmann, K., “Force and Pressure Measurements 
Including Surface Row Visualization on a Cropped 
Delta Wing” Proceedings of Symposium on Interna- 
tional Vortex Flow Experiment on Euler Code Vali- 
dation, Stockholm, Sweden, October 1-3, 1986, pp. 
63-87. 

25. Butefisch, K. A., Pallek, D. and Sauerland, K., - 
H., “International Vortex Row Experiment-Results of 
Three Component LDA Measurements on a 65° Delta 
Wing.” DFVLR IB 222-87 A 34, 1987. 

26. Elsenaar, A., Hjelmberg, L.. Butefisch, K. and Ban- 
ning W. J., "The International Vortex Row Experi- 
ment,” AGARD-CP-437, Lison, Portugal, May 1988, 
VoL 1., pp. 9.1-9.23. 

27. Hitzel, S. M., “Wing Vortex-Rows Up into Vortex- 
Breakdown-A Numerical Simulation," AIAA 88- 
2518-CP, 1988, pp. 73-83. 

28. Bannink. W. J. and Houtman, E. M., “Experimental 
and Computational Study of the Vortical Row Over 
a Delta wing at High Angles of Attack,” IUTAM 
Symposium on Ruid Dynamics of High Angle of 
Attack, University of Japan, Tokyo, Japan, September 
14-17, 1992. 

29. Laine, S., Siikonen, T. and Kaurinkoski, P-, “Cal- 
culation of Transonic Viscous Row Around a Delta 
Wing,” ICAS 92-4.2.1, Beijing, PRC. September 22- 
25. 1992, pp. 286-295. 

30. Erickson, G. E„ “Wind Tunnel Investigation of The 
Interaction and Breakdown Characteristics of Slender- 
Wtng Vortices at Subsonic, Transonic and Supersonic 
Speeds,” NASA Tech, paper 3114, November 1991. 


5 



sm 





ST 























Sp*n-« t •• dl»t*ne». 



Fig. 5. 


Static-pressure and Mach contours on the wing and the 
plane of symmetry; Moo = 0.85, a = 20°, t = 3.6 
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Fig. 9. Static-pressure contours on the wing and the plane of 
symmetry; M» = 0.85, a = 20°, t = 5.52 
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Fig. 10. Ray lines on the wing surface and the static-pressure 
variation along them; = 0.85, a = 20°, t = 5.52 
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Fig. 11. Total-Mach contours, streamlines and velocity vectors on a ray plane passing 
through the vortex-breakdown bubbles; Moo = 0.85, a = 20°, t = 5.52 
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Fig. 12 Surface-pressure and Mach contours and particle trance on wing 
and symmetry planes; M 0 0 = 0.85, a = 20°, t = 3.6 
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Supersonic Vortex Breakdown on a Delta Wing 
M = 0.85, Re - 3,230,000 and AOA = 20 
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Fig. 13. Surface-pressure contours and total-pressure surfaces for a wing plan view; Mach contours 
and total-pressure surface on a ray plane passing through the vortex-breakdown bubble; 
M x = 0.S5. q = 20°. t = 3.6 


13 



tbf 


AIAA 93-0526 

THREE-DIMENSIONAL SUPERSONIC 
VORTEX BREAKDOWN 


Osama A. Kandil and Hamdy A. Kandil 
Old Dominion University, Norfolk, VA 23529 

C. H. Liu 

NASA Langley Research Center, Hampton, VA 23665 


31st Aerospace Sciences 
Meeting & Exhibit 

January 11-14, 1993 /Reno, NV 


For permission to copy or republish, contact the American Institute of Aeronautics and Astronautics 
370 L'Enfant Promenade, S.W., Washington, D.C. 20024 





THREE-DIMENSIONAL SUPERSONIC VORTEX BREAKDOWN 

Osama A. Kandil* and Hamdy A. Kandil** 

Old Dominion University, Norfolk, VA 23529 
and 

C. H. Liu*** 

NASA Langley Research Center, Hampton, VA 23665 


ABSTRACT 

Three-dimensional, supersonic vortex -breakdown prob- 
lems in bound and unbound domains are solved. The 
solutions are obtained using the time-accurate integra- 
tion of the unsteady, compressible, full Navier-Stokes 
(NS) equations. The computational scheme is an implicit, 
upwind, flux-difference splitting, finite-volume scheme. 
Two vortex-breakdown applications are considered in the 
present paper. The first is for a supersonic swirling jet 
which is issued from a nozzle into a supersonic uniform 
flow at a lower Mach number than that of the swirling jet 
The second is for a supersonic swirling flow in a config- 
ured circular duct. In the first application, an extensive 
study of the effects of grid fineness, shape and grid-point 
distribution on the vortex breakdown is presented. Four 
grids are used in this study and they show a substantial 
dependence of the breakdown bubble and shock wave on 
the grid used. In the second application, the bubble-type 
and helix-type vortex breakdown have been captured. 

INTRODUCTION 

Longitudinal vortex/transverse shock-wave interac- 
tion is a complex flow phenomenon which develops in 
several external and internal flow applications. For exter- 
nal flows, the transonic flow around a delta wing in the 
high-angle-of-attack range 1,2 and the transonic and super- 
sonic flows around a strake-delta wing configuration in 
the moderate to high-angle-of-attack range 3 are some of 
the applications. Vortex-breakdown usually occurs behind 
the transverse shock wave over the delta wing resulting in 
a loss of lift. Such a breakdown is undesirable and flow- 
control methods need to be developed to eliminate the 
vortex breakdown. For internal flows, the supersonic in- 
let ingesting a vortex and the supersonic combustor where 
fuel is injected in a swirling jet to enhance fuel-air mixing 
are some of the applications. Vortex breakdown behind 
the transverse shock wave in these applications is desir- 
able since it enhances mixing and stability of the flame 4 " 7 , 
and hence its occurrence need to be controlled for opti- 
mum performance. 
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For such problems, computational schemes are needed 
to study, predict and control vortex-shock interaction in- 
cluding vortex breakdown. Unfortunately, the literature 
lacks this type of analysis with the exception of the pre- 
liminary work of Liu, Krause and Menne 8 , Copening and 
Anderson 3 , Delery, ei al. 4 and Kandil and Kandil 9 . 

The first time-accurate NS solution for a quasi- 
axisymmetric supersonic vortex-breakdown was devel- 
oped by Kandil, Kandil and Liu in Ref. 10. A supersonic 
quasi-axisymmetric vortex flow in a configured circular 
duct was considered. The time-accurate solution of the 
unsteady, compressible NS equations was obtained us- 
ing an implicit, upwind, flux -difference splitting finite- 
volume scheme. A shock wave was generated near the 
duct inlet and unsteady vortex -breakdown was predicted 
behind the shock. The predicted flow was character- 
ized by the evolution, convection and shedding of vortex- 
breakdown bubbles. The Euler equations were also used 
to solve the same problem. The Euler solution showed 
larger sfce and number of vortex-breakdown bubbles in 
comparison with those of the NS solutions. The time- 
accurate solution was carried out for 3,200 time steps 
which were equivalent to a dimensionless time of 16. 
Only one value of Reynolds number of 10,000 was con- 
sidered in Ref. 10. 

In a later paper by Kandil, Kandil and Liu 11 , the 
study of this flow was extended using time-accurate com- 
putations of the NS equations with a fine grid in the 
shock-vortex interaction region and for long computa- 
tional times. Several issues were addressed in that study. 
First, the effect of Reynolds number on the temporal evo- 
lution and persistence of vortex-breakdown bubbles be- 
hind the shock was shown. In that stage of computations, 
the conditions at the downstream exit were obtained by 
extrapolating the components of the flowfield vector from 
the interior cell centers. Although the flow was super- 
sonic over a large portion of the duct exit, subsonic flow 
existed over a small portion of the exit around the duct 
centerline. Therefore, selected flow cases were computed 
using Rieraann-invariant-type boundary conditions at sub- 
sonic points of the duct exit. The effect of swirl ratio at 
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the duct inlet was also investigated. Recently, the critical 
effects of the downstream boundary conditions on the su- 
personic vortex -breakdown was extensively investigated 
by the same authors 12 for both internal and external flows. 

In Refs. 10-12, the present authors assumed the flow 
quasi-axisymmetric and the NS equations were solved 
using the three-dimensional solver "FTNS3D" by forc- 
ing the components of the flowfield vector to be equal 
on two axial planes in close proximity of each other. 
Quasi-axisymmetric solutions are one-order of magnitude 
less in computational cost than the corresponding three- 
dimensional solutions, and they still provide substantial 
physical understanding of the supersonic vortex break- 
down. At this substantially reduced computational cost, 
we were able to study the effects of Reynolds number, 
swirl ratio and downstream exit-boundary conditions us- 
ing time-accurate stepping. However, the previous exper- 
imental studies for both incompressible 1 ** 15 and super- 
sonic vortex breakdown 16 showed that the flow is three 
dimensional. 

Hence, we consider the three-dimensional solution 
of the NS equations to realistically simulate the vor- 
tex breakdown problem. The computational solution of 
two main vortex-breakdown problems are presented in 
this paper the first is for vortex breakdown of a super- 
sonic swirling jet issued from a nozzle into a supersonic 
uniform flow at a lower Mach number than that of the 
swirling jet, and the second is for vortex breakdown of a 
supersonic swirling flow in a configured circular duct In 
the first problem, an extensive study of the effects of grid 
fineness, shape and grid-point distribution on the break- 
down bubble is presented. Four grids have been used in 
the study, and they show a substantial dependence of the 
breakdown bubble and shock wave on the grid used. In 
the second problem, the time-accurate solution shows two 
modes of supersonic vortex breakdown; a bubble type 
and a spiral type. 

HIGHLIGHTS OF THE FORMULATION 
AND COMPUTATIONAL SCHEME 

The conservative, unsteady, compressible, lami- 
nar full Navier-Stokes equations in terms of time- 
independent, body-conformed coordinates , £ 2 and 
are used to solve the problem. The equations are given 
in Ref. 11 and hence they are not shown here. Along 
with these equations, boundary conditions are specified 
at the computational-domain inlet, side wall and down- 
stream exit. The boundary conditions are presented in the 
next section. The initial conditions will also be presented 
in the next section. 

The computational scheme used to solve the unsteady, 
compressible full NS equations is an implicit, upwind, 
flux-difference splitting, finite-volume scheme. It em- 
ploys the flux-difference splitting scheme of Roe which is 
based on the solution of the approximate one-dimensional 
Riemann problem in each of the three directions. In the 


Roe scheme, the in viscid flux difference at the interface 
of a computational cell is split into left and right flux dif- 
ferences. The splitting is accomplished according to the 
signs of the eigenvalues of the Roe averaged-Jacobian 
matrix of the in viscid flux at the cell interface. The 
smooth limiter is used to eliminate oscillations in the 
shock region. The viscous and heat-flux terms are lin- 
earized and the cross-derivative terms of the viscous Ja- 
cobians are dropped in the implicit operator. These terms 
are differenced using second-order spatially accurate cen- 
tral differencing. The resulting difference equation is ap- 
proximately factored and is solved in three sweeps in the 
£>, f 2 and ( 3 directions. The scheme is coded in the 
computer program which is called **FTNS3D M . 

COMPUTATIONAL RESULTS AND DISCUSSION 

I. Three-Dimensional Vortex Breakdown of 
a Supersonic Swirling Jet 

A supersonic swirling jet at a Mach number of M, = 

3.0 is issued from a nozzle into a supersonic uniform flow 
at a Mach number of M x = 2.0. The freestream Reynolds 
number, R „ is 296,000. The nozzle-exit radius is the 
characteristic length and the length of the computational 
domain is 7.0 dimensionless units. The purpose of the 
present computational case is to simulate the flow of the 
experimental study of Ref. 16. It was reported in Ref. 16 
by Metwally, et al. that it was difficult to detect any 
vortex-breakdown bubble behind the formed shock in the 
swirling jet flow and that the shock was oscillating around 
a mean position. For the present computational study, it 
is decided to use four types of structured grids to solve the 
unsteady, compressible NS equations accurately in time. 
The cross-section of the computational domain is taken as 
a square section for three types of grid which are called 
Grid type 1, 2 and 3. For the fourth grid, Grid type 4, a 
circular section is used. For Grid types 1, 2 and 3, the 
length of one-half the square-section side is 3.5 units and 
for Grid type 4, the radius length of the circular section 
is 3.5 units. A time step of 0.02 is used for all the four 
types of grids. 

1.1 Boundary and Initial Conditions: The inflow bound- 
ary conditions are adapted from the limited experimen- 
tal data of Ref. 16. Unfortunately, the experimental data 
available in Ref. 16 are given along one diameter only 
of the circular section. The profiles of the experimen- 
tal data are not symmetric with respect to the diameter 
center point To produce three-dimensional profiles from 
the experimental data, two methods are used. In the first 
method, the asymmetry of the experimental profiles is 
maintained by assuming the profiles on the right-hand 
side of the initial cross-section to be the same as those of 
the upper half of the experimental data and the profiles 
on the left-hand side of the initial cross-section to be the 
same as those of the lower half of the experimental data. 

In the second method, the initial cross-section profiles are 
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assumed to be quasi-axi symmetric and are set equal to the 
profiles of the upper half of the experimental data. The 
average swirl ratio of the asymmetric profiles is 0.2 and 
the swirl ratio of the quasi-axi symmetric profiles is 0.22. 
The experimental data are used for r < 1, and for 1 < 
r S 3.5 uniform wind-tunnel conditions corresponding to 
Moc = 2.0 are used. Figures 1 and 2 show the generated 
asymmetric and quasi-axisymmetric distributions of the 
axial velocity, u, tangential velocity v, density, p, and 
pressure, p. 

The boundary conditions on the outer boundaries of 
the computational domain are assumed to be uniform 
conditions corresponding to M x = 2.0. 

The outflow boundary conditions at the exit of the 
computational domain are obtained by extrapolating all 
the flowfield components from the interior cell centers 
at the exit This is justified since the experimental data 
of Ref. 16 show that the flow becomes supersonic a few 
radii downstream of the shock/vortex interaction region. 

The initial conditions in the entire domain correspond 
to frees tream conditions at - 2.0. Hence, the present 

flow case simulates a sudden discharge of a swirling 
supersonic jet from a nozzle into a uniform supersonic 
flow. 

1.2 Grid Type 1: A rectangular grid of 210x51x51 

points in the axial (x direction) and cross-flow directions 
(y and z directions), respectively is considered. The grid 
points are clustered in the axial direction near the inflow 
boundary and near the vortex-core axis. The minimum 
grid size is 0.057 in the y and z directions and 0.0147 
in the x direction. Figure 3 shows a side view and a 
cross-flow plane of the grid. 

In Fig. 4, snapshots of the streamlines on a horizon- 
tal plane passing through the computational domain axis 
and the total Mach contours on the same plane are given 
at selected dimensionless times. At t = 2.0, the Mach 
contours show the existence of a strong shock wave at 
the centerline. The shock is formed due to the mismatch 
between the static pressures of the supersonic swirling 
jet and the supersonic uniform surrounding flow. Behind 
the shock wave, the streamlines indicate the existence 
of a reversed flow region. At t = 3.0, the recirculating 
bubble-flow grows in size and moves upstream. The for- 
mation of a two-bubble cell is clear. The Mach contours 
show that the shock also moves upstream ahead of the 
two-bubble cell. At t = 4.0 and 5.0, the shock wave and 
the two-bubble cell still move upstream while the bubbles 
are growing in size. The bubbles reach their maximum 
size at t = 5.0. The snapshots of the solution for t £ 10.0 
show that the bubble-shock system is quasi-steady and 
oscillating around a mean axial location. The maximum 
change in the bubble size is less than 10%. The axial ve- 
locity recovers its supersonic values at the exit boundary, 
and hence the use of extrapolation boundary conditions 
is justified. 


L3 Grid Type 2: This grid consists of 145x61x61 

points in the axial and cross-flow directions, respectively. 
The grid points are redistributed to obtain better solution 
of the vortex core and the vortex-shock interaction region. 
The minimum grid size is 0.024 in the y and z directions 
and 0.014 in the axial direction. Figure 5 shows a side 
view and a cross-flow plane of the grid. 

In Fig. 6, snapshots of the streamlines and total Mach 
contours on a horizontal plane passing through the com- 
putational domain axis are given at selected dimension- 
less times. At l = 2.0, a two-bubble cell is formed behind 
a sharply oblique conical shock. The bubbles sizes are 
larger than those of Grid type 1 at the same time level. 
At t = 3.0, the two-bubble cell grows in time and moves 
upstream along with the shock wave. Again, the bubbles 
are larger in size and closer to the inflow boundary than 
those of Grid type 1. Moreover, it is interesting to notice 
that the bubbles are longer in the axial direction than their 
length in the lateral direction in comparison with those of 
Grid type 1. It should be noted here that the number 
of grid points around the computational domain axis for 
the present grid is larger than those of Grid type 1 and 
the number of grid points in the axial direction is less 
than those of Grid type 1. At t = 4.0, the two-bubble cell 
moves downstream, another small bubble appears and the 
shock splits into two shocks; a weak shock which is fol- 
lowed by a strong shock surrounding the bubbles. At t 
= 5.0, the shock shape changes to become more oblique 
and the two-bubble cell grows slightly. As the solution is 
advanced in time, it is observed that continuous changes 
in the bubbles size, shape and location occur with larger 
amplitudes than those of Grid type 1. For t S 10.0, the 
bubbles show highly unsteady flow with an oscillating 
shock wave. 

1.4 Grid T>pe 3: In this grid, the number of grid points 
is kept the same as that of Grid type 2. The grid points are 
redistributed in the axial and cross-flow directions. In the 
axial directions, 90 grid points are used in the range of x 
= 0 to 2.0, in comparison with 71 grid points in Grid type 
2. The minimum cell size in the axial direction is 0.0084. 
In the cross-flow plane, the grid points are redistributed 
such that the grid aspect ratio does not exceed 4.0. Figure 
7 shows the grid and Fig. 8 shows the results. 

The results show that at t = 2.0, a small bubble is 
captured off the computational domain axis. Later on, 
at t = 4.0, the small bubble disappears. Another bubble 
is captured at t = 6.0 and it also disappears later on. 
For t > 8.0, no more bubbles are captured. A strong, 
al most-normal shock is captured around the axis. It is 
located a little more downstream from the inflow section 
in comparison with that of Grid type 2. 

Next, it is decided to use the quasi-axisymmetric ini- 
tial profiles (Fig. 2) which have higher swirl ratio than 
those of the asymmetric initial profiles (Fig. 1). The re- 
sults are shown in Fig. 9. The results show the formation 
of a small two-bubble cell. The bubbles shape changes 
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slowly with time and the shock oscillates with very small 
am plitude, around a mean position which is a little more 
downstream than that of the previous case and the case 
of Grid type 2. 

L5 Grid Type 4: A circular grid consisting of 145 x61 x 
49 points in the axial, radial and circular directions, re- 
spectively, is used. The grid points are clustered around 
the axis for a good resolution of the vortex core and 
around r * 1 for good resolution of the shear layer be- 
tween the swirling jet and the uniform freestream flow. 
In the axial direction, the grid points are distributed as 
those of Grid type 3. The circular grid has the advan- 
tage of offering better resolution near the axis, where it 
is needed. Moreover, with the circular grids, the number 
of grid points along a diameter in the cross-flow plane 
is doubled without increasing the total number of grid 
points, in comparison with the previous grids of square- 
section cross-flow planes. Figure 10 shows the circular 
grid. 

As in the case of Grid type 3, two sets of initial 
profiles; namely the quasi-axisymmetric and asymmetric 
profiles, are used with this grid. As with Grid type 3 
and for the asymmetric initial profiles, a small bubble is 
formed behind the shock and disappears after a few time 
steps. 

Figure 11 shows snapshots of the streamlines and 
Mach number contours for the quasi-axisymmetric initial 
profiles. The results show the formation of ajnulti- 
bubble vortex breakdown behind the central strong part of 
the shock system. A two-bubble cell is then established 
and persists for the rest of the computational time. The 
relative size of the two bubbles is continuously changing 
and the global picture is looked at as a quasi-steady one. 

This study exclusively shows why it was very difficult 
to see any vortex-breakdown bubbles as was reported in 
Ref. 16. It is understood now in view of the results of 
the four grids that the size of the bubbles are either very 
small to be seen for quasi-axisymmetric initial profiles or 
they are transient bubbles for asymmetric initial profiles. 

II. Three-Dimensional Supersonic Vortex 

Breakdown in a Configured Duct 

The computational domain consists of a configured 
circular duct with a total length of 2.9 dimensionless units, 
where the duct radius is used as a characteristic length. 
The duct consists of a constant diameter cylindrical por- 
tion of radius one followed by a divergent portion which 
is intended to stabilize the formed shock wave, a constant 
cylindrical part and finally a convergent-divergent nozzle 
which is intended to accelerate the exhaust flow to su- 
personic speeds. The grid consists of 200x51x49 grid 
points in the axial, radial and circumferential directions, 
respectively. The grid points are clustered near the inlet 
section in the axial direction for a good resolution of the 
shock and the shock/vortex interaction region, and in the 


cross-flow plane around the duct axis for a good resolu- 
tion of the vortex core. The gird points are also clustered 
near the duct walls for a good resolution of the bound- 
ary layer. The minimum cell size is 0.002. Figure 12 
shows the computational grid. The freestream conditions 
correspond to M x ■ 1.75 and Reynolds number of 10 s . 

n.l Boundary and Initial Conditions: The initial pro- 
file for the tangential velocity is given by 



where Uoc = 1.74, r m = 0.2 and k, = 0.1. The max- 
imum j?-, swirl ratio /?, is at r = 0.224 and its value is 
kept at 032. The radial velocity, w, at the initial station 
is set equal to zero and the radial momentum equation is 
integrated to obtain the initial pressure profile. Finally, 
the density p is obtained from die definition of the speed 
of sound for the inlet flow. With these compatible set 
of profiles, the computations are carried out accurately in 
time with At = 0.0025 for two computational applica- 
tions. The duct-wall boundary conditions follow the typ- 
ical Navier-Stokes solid-boundary conditions for the first 
case. For the second case, in viscid duct-wall boundary 
conditions are used to reduce the effect of the boundary- 
layer separation on the vortex breakdown process. The 
downstream exit-boundary conditions are obtained by ex- 
trapolating the flowfield components from the interior cell 
centers at the exit. 

The initial conditions correspond to stagnation condi- 
tions throughout the interior computational domain. 

D.2 Viscous Duct-Wall Boundary Conditions: Figure 
13 shows snapshots of the solution at selected time levels. 
At t = 2.0, a small recirculating region is formed behind 
the strong normal part of the shock wave. Two stagna- 
tion points are recognized along the duct axis. The total 
Mach number contours show the position of the shock 
front near the inflow section and the position of the re- 
circulation zone behind the shock wave. As the com- 
putations is advanced in time, the bubble size grows in 
the axial and radial directions and the shock-bubble sys- 
tem moves downstream. At t = 3.5, it is noticed that the 
bubble size grows and the shock wave is deformed ac- 
cordingly. The solution is quasi-axisymmetric as shown 
by the streamlines and Mach number contours. Starting 
at t = 4.0, the bubble size grows in the lateral direction, 
moves upstream towards the inflow boundary pushing the 
shock wave upstream. Small flow asymmetry is also no- 
ticed. For t 2 5.5, another phase of the solution history 
develops, where a reversed normal shock is formed inside 
the vortex-breakdown bubble (see Mach contours at t ■ 
6.0). The normal shock wave turns the reversed flow to 
subsonic. As the computations advance in time, the bub- 
ble system starts to move downstream towards the duct 
exit with a new recirculating region formed behind the 
shock wave. The flow becomes quite asymmetric. 
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At t = 11.5, as the computations advance in time, 
some of the features noticed experimentally for incom- 
pressible vortex flows in pipes could be recognized; e.g., 
the asymmetric vortex breakdown seen at t = 11.5. The 
streamlines clearly show the spiral-type of vortex break- 
down and the asymmetric shedding of the vortex break- 
down bubbles. It should be noted here that such phe- 
nomenon was not captured using the quasi-axisymmetric 
assumption in Ref. 11 by the present authors. The shed- 
ding of the vortex-breakdown bubbles continues as new 
bubble systems were formed behind the shock wave. At 
t = 16.5, the shedding of two asymmetric bubbles can be 
seen and a two-bubble system is formed upstream along 
with a very small recirculating region just downstream 
of the central part of the shock wave. At t = 19, it is 
observed that two bubbles merge together while they are 
convected and shed. 

An important parameter affecting the flow in the duct 
is the interaction of the shock with the boundary-layer 
flow on the duct wall. This causes the separation of 
the wall boundary layer, as can be seen at t = 31, 42 
and 43. The Mach number contours show the separation 
of the boundary as a result of the interaction with the 
shock wave, and the streamlines show the reduction in the 
vortex-breakdown-bubble size as a result of the boundary- 
layer thickening. The shedding of the inclined vortex 
rings shown at t = 33 is similar to the spiral type of 
vortex breakdown, where the upper part of the vortex 
rings rotates in the clockwise direction and the lower part 
rotates in the opposite direction. A new vortex ring is 
formed behind the shock while the spiral-like system is 
moving downstream. 

As the solution is advanced in time, the size of the 
breakdown region is reduced in the radial direction as 
a result of the boundary-layer thickening. An interesting 
breakdown mode is shown at t = 38.5 where a bubble-type 
vortex breakdown followed by a spiral-type is formed 
downstream of the shock wave. This phenomenon was 
observed in the experimental studies of Sarpkaya 13 and 
has never been captured computationally. 

The reduction of the breakdown-region size continues 
with the advance in the time as can be seen at t = 41, 42 
and 44. At t = 46, no recirculation zone is observed 
and the vortex-breakdown system is dissipated totally. 
It is also observed that the shock-wave system moves 
continuously in the downstream direction as of t > 31. 

11.3 Inviscid Duct-Wall Boundary Conditions: The 

effect of the shock/boundary-layer interaction at the duct 
wall is further investigated by treating the duct wall as an 
inviscid wall. At l = 43.0, inviscid- wall boundary condi- 
tions are applied at the duct wall with ail the other bound- 
ary conditions remaining the same. Samples of the results 
are shown in Fig. 14. At t = 43.5, the vortex-breakdown 
bubbles are recovered and the shock wave becomes nor- 
mal to the duct wall. The shedding of the vortex rings 
continues as the solution is advanced in time as can be 


seen at t = 45.5, where the vortex rings are recognized. 

It is also noticed that the vortex-breakdown bubble size 
starts to increase in the radial direction. Further increase 
in the breakdown-region size is noticed at t ■ 47. It is also 
noticed that the position of the shock wave with respect to 
the duct inlet is fixed while the shape of the central part 
is continuously changing according to the shape of the 
bubbles behind the shock. The shedding of the vortex- 
breakdown bubbles continues in an asymmetric form as 
can be seen at t = 63, 69 and 72 5 . It is interesting to 
notice that the vortex-breakdown system survives and is 
not dissipated as in the case of viscous duct wall. The 
computations is advanced until t = 75 without any sign 
of dissipation of the vortex-breakdown-bubbles. It is con- 
cluded that the disturbances caused by the wall boundary- 
layer separation are the reason behind the disappearance 
of the vortex-breakdown system in the case of viscous 
duct wall. This might be caused by the pressure gradi- 
ents resulting from the change in the vortex-core outer 
boundaries. 

Sarpkaya 13 noticed that the boundary-layer separation 
and reversed flow occurred on the tube wall in the case 
of a swirling incompressible flow in a divergent tube. He 
suggested that the bubble pressure gradient caused by the 
tube divergence and that caused by the vortex breakdown 
are the reasons behind the separation. He concluded that 
the viscous effects on vortex breakdown in tubes are very 
significant and that, because of the flow separation, a bet- 
ter simulation of the vortex breakdown is not likely to 
emerge from solving numerically the full Navier-Stokes 
equations even if the problems of numerical instability 
were to be solved. In the case of supersonic vortex break- 
down, the problem is much more involved because of the 
shock/boundary layer interaction and the assumption of 
inviscid walls seems to isolate the wall viscous effects. 

CONCLUDING REMARKS 

Three-dimensional, supersonic vortex-breakdown flows 
in bound and unbound domains are simulated compu- 
tationally using the time-accurate solution of the un- 
steady, compressible, laminar, full Navier-Stokes equa- 
tions. Two main vortex-breakdown applications are con- 
sidered in this paper. The first application is for a su- 
personic swirling jet issued in a supersonic uniform flow 
at a lower Mach number. This flow case was consid- 
ered earlier by Metwally and his co-workers in Ref. 16, 
where it was reported that no vortex breakdown bubble 
was seen behind the shock-wave system. A systematic 
computational investigation was carried out using four 
types of grids which ranged from coarse- to find-grid dis- 
tributions and from rectangular to circular grid lines in 
the cross-flow planes. It has been shown that the coarse 
grid produces large vortex bubbles and the fine girds pro- 
duce either transient small vortex bubbles or quasi-steady 
small vortex bubbles. Using the fine-grid results as the 
ones closely representing the experimental flow case, it 
is understood why Metwally, et fid. 16 were not able to see 
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any vortex-breakdown bubble. The bubbles were either 
small or small and transient Moreover, this study shows 
why previous investigators 17,11 of incompressible vortex- 
breakdown flows w ere able to produce numerical results 
using coarse grids (coa r ser than die coarse grid used in the 
present paper) at low Reynolds numbers which were com- 
parable to experimental results at high Reynolds numbers. 
It is now understood that coarse grids have made it pos- 
sible to simulate experimental results at high Reynolds 
number. 

The second application is for a supersonic swirling 
flow in a configured circular duct Here, the duct-wall 
boundary conditions are used once for a viscous wall 
and another for an in vise id wall. With the viscous-wall 
boundary conditions, it has been observed that the vor- 
tex breakdown is transient and it has been dissipated by 
the effect of separated flow from the duct-wall bound- 
ary layer. However, during the transient formation of 
vortex-breakdown flows, both the bubble-type and spiral- 
type vortex breakdown are captured. Spiral-type vortex 
breakdown was not captured in Ref. 11 by the present 
authors due to the quasi-axisymmetric assumption used. 
With the inviscid duct-wall boundary conditions, the vor- 
tex breakdown is persistent and does not dissipate. The 
three-dimensional relieving effect on the vortex break- 
down modes is apparent from the present results when 
they are compared with those of Ref. 11, where the quasi- 
axisymmetric assumption has been used. 
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Fie. 5. Grid type 2 (rectangular fine gird in the cross-flow plane), 

145x61 x61 grid points in the axial and cross-flow plane, respectively. 
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Fig. 7. Grid type 3 (rectangular fine grid), 145x61x61 grid points in the 
axial and cross-flow plane, respectively. 
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Fig. 10. Grid* type 4 (circular fine grid), 145x61x61 grid points in the 
axial and cross-flow plane, respectively. 
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Fig. 1 1 . Streamlines and Mach contours in a horizontal plane for a super- 
sonic swirling jet using Grid type 4. M, = 3.0, M* = 2.0 and R* = 296,000. 
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Fig. 13. Streamlines and Mach contours in a horizontal plane for a super- 
sonic swirling flow in a circular duct, = 1.75, 0 = 0.32 and 
R, = 100,000, viscous wall. 
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Fig. 14. Streamlines and Mach contours in a horizontal plane for a super- 
sonic swirling flow in a circular duct, M Q 0 = 1.75, 0 = 0.32 and 
R< ~ 100,000, inviscid wall. 
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Abstract — Existing numerical simulations and physical aspects of subsonic and supersonic vortex-break- 
down modes are reviewed. The solution to the problem of supersonic vortex breakdown is emphasized 
in this paper and carried out with the full Navier-Stokes equations for compressible flows. Numerical 
simulations of vortex-breakdown modes are presented in bounded and unbounded domains. The effects 
of different types of downstream-exit boundary conditions are studied and discussed. 


INTRODUCTION 

At sufficiently high angles of attack, vortex breakdown has been observed along the primary 
leading-edge vortex cores of a delta wing. Two distinct forms of vortex breakdown have been 
documented experimentally [1]. The first form is the bubble type, and the second is the spiral type. 
The bubble type shows a sudden, nearly axisymmetric swelling of the vortex core into a bubble; 
the spiral type shows an asymmetric spiral vortex filament, which is followed by a rapidly spreading 
turbulent flow. Both types of breakdown are characterized by an axial stagnation point and a 
limited region of reversed axial flow. Much of our knowledge of vortex breakdown has been 
obtained from experimental studies of incompressible flows in pipes, where both types of 
breakdown (and other types as well) are generated and documented [2-4], 

The major effort in the numerical studies of vortex-breakdown flows has been focused on 
incompressible, quasi-axisymmetric isolated vortices. Grabowski and Berger [5] were the first to use 
the incompressible, quasi-axisymmetric Navier-Stokes (NS) equations to study isolated vortex flow 
in an unbounded region. Hafez et al. [6] solved the incompressible, steady, quasi-axisymmetric 
Euler and NS equations with the stream function-vorticity formulation for isolated vortex flows. 
They predicted vortex-breakdown flows similar to those of Garbowski and Berger [5]. Salas and 
Kuruvila [7] solved the unsteady, quasi-axisymmetric NS equations in a straight circular pipe and 
obtained steady, multiple-bubble vortex breakdown for the Reynolds number, Re (based on the 
pipe diameter), range 100-1800. Menne [8] has also used the stream function-vorticity formulation 
for studying unsteady, incompressible quasi-axisymmetric isolated vortex flows. Wu and Hwang [9] 
used the stream function-vorticity formulation to study quasi-axisymmetric vortex breakdown in 
a pipe. Their study focused on the effects of inflow conditions, wall boundary conditions and Re 
on breakdown structure. They showed that the evolution of breakdown can be steady, periodic 
or unsteady, dependent on the inflow velocity profiles and Re. Menne and Liu [10] integrated the 
laminar, incompressible NS equations for breakdown of a vortex in a slightly diverging pipe. They 
showed breakdown flow cases that are based on the purely quasi-axisymmetric and nonaxisymmet- 
ric analyses. The results were in good agreement with the experimental results of Leibovich [4]. Spall 
et al. [11] used the vorticity-velocity formulation of the incompressible NS equations to study 
three-dimensional vortex breakdown. Breuer and Hanel [12] solved the unsteady incompressible 
NS equations with a dual time-stepping, upwind scheme to study the temporal evolution of the 
three-dimensional vortex breakdown. In Refs [11, 12], both types of breakdown (the bubble and 
the spiral type) were predicted. Reviews of the physical and computational aspects of the 
incompressible vortex breakdown have been presented by Krause [13, 14]. One of the most 
important aspects of vortex breakdown, which Krause discusses in Ref. [14], is the effect of side-wall 
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boundary conditions on the upstream or downstream motion of the breakdown point. Also, he 
discusses different outflow boundary conditions. 

Longitudinal vortex and transverse shock-wave interactions are typical phenomena that appear 
in transonic and supersonic flows over strake-wing configurations at moderate to high angles of 
attack, at a supersonic inlet as a vortex is injested and inside a supersonic combustor where fuel 
is injected in a swirling jet to enhance fuel-air mixing [15-17], For the strake-wing configuration, 
vortex breakdown is undesirable because the stall phenomenon results; hence, its occurrence must 
be delayed. On the other hand, vortex breakdown for the other two applications is desirable 
because it enhances the mixing of air and fuel and stability of the flame [18, 19]; hence, its 
occurrence must be controlled for optimum performance. Unfortunately, the literature lacks this 
type of analysis, with the exception of the preliminary work of Liu et al. [29], Copening and 
Anderson [21], Delery et al. [15], Kandil and Kandil [22] and Meadows et al. [23], 

The first time-accurate, NS solution for a supersonic vortex breakdown was developed by Kandil 
et al. [24], A supersonic quasi-axisymmetric vortex flow in a configured circular duct was 
considered. The time-accurate solution to the unsteady, compressible NS equations was obtained 
with an implicit, upwind, flux-difference splitting, finite-volume scheme. A shock wave was 
generated near the duct inlet and unsteady vortex breakdown was predicted behind the shock. The 
predicted flow was characterized by the evolution, convection and shedding of vortex-breakdown 
bubbles. The Euler equations were also used to solve the same problem. The Euler solution 
predicted increases in both the size and number of vortex-breakdown bubbles, in comparison with 
the NS solutions. The time-accurate solution was carried out for 3200 times steps, which was 
equivalent to a dimensionless time of 16. A single Re value (10,000 based on the inlet radius) was 
considered in that study [24], 

In a later paper [25], the study of this flow was extended with time-accurate computations of 
the NS equations with a fine grid in the shock-vortex interaction region and for longer 
computational times. Several issues were addressed in that study. First, the effect of Re on the 
temporal evolution and persistence of vortex-breakdown bubbles behind the shock was shown. In 
that stage of the computations, the conditions at the downstream exit were obtained by 
extrapolating the components of the flow-field vector from the interior cell centers. Although the 
flow was supersonic over a large portion of the duct exit, subsonic flow existed over a small portion 
of the exit around the duct centerline (CL). Therefore, selected flow cases were computed with 
Riemann-invariant boundary conditions at the subsonic points of the duct exit. The effect of swirl 
ratio at the duct inlet was also investigated. 

Recently, the critical effects of the downstream boundary conditions on the supersonic vortex 
breakdown have been investigated extensively by the same authors [26] for both internal and 
external flows. 

In the present paper, the numerical simulation of supersonic vortex-breakdown flows for 
bounded and unbounded domains are reviewed. The effects of Re, swirl ratio and downstream-exit 
boundary conditions are considered and discussed along with certain physical and numerical issues. 


OVERVIEW OF THE FORMULATION AND COMPUTATIONAL SCHEME 

The conservative, unsteady, compressible, full NS equations, in terms of the time-independent, 
body-conformed coordinates £ 2 and have been used to solve the problem. The equations are 
given in Ref. [25] and are not shown here. With these equations, boundary conditions are specified 
at the computational domain inlet, side wall and downstream exit. The downstream-exit boundary 
conditions are presented and discussed in the next section. The initial conditions will also be 
presented in the next section. 

The computational scheme used to solve the unsteady, compressible, full NS equations is an 
implicit, upwind, flux-difference splitting, finite-volume scheme. This scheme employs the flux- 
difference splitting scheme of Roe, which is based on the solution to the approximate one-dimen- 
sional Riemann problem in each of the three directions. In the Roe scheme, the inviscid flux 
difference at the interface of a computational cell is split into left and right flux differences. The 
splitting is accomplished in accordance with the signs of the eigenvalues of the Roe averaged- 
Jacobian matrix of the inviscid flux at the cell interface. The smooth limiter is used to eliminate 
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oscillations in the shock region. The viscous and heat-flux terms are linearized and the cross- 
derivative terms of the viscous Jacobians are dropped in the implicit operator. These terms are 
differenced with second-order, spatially-accurate central differencing. The resulting difference 
equation is approximately factored and is solved in three sweeps in the <* 2 and f 3 directions. 
The scheme is coded in the computer program “FTNS3D”. 

The quasi-axisymmetric solutions are obtained with the three-dimensional code by forcing the 
flow-field vectors to be equal on two axial planes, which are in close proximity to each other. 
Quasi-axisymmetric solutions require 1 order of magnitude less in computational time than the 
three-dimensional solutions. They still provide substantial physical understanding of the supersonic 
vortex breakdown and the dominant parameters that affect it. 


COMPUTATIONAL RESULTS AND DISCUSSION 


Vortex Breakdown in a Configured Circular Duct 

Figure 1 shows a configured circular duct with a short, straight cylindrical part at the inlet which 
is followed by a short divergent cylindrical part until the axial length of 0.74 dimensionless unit, 
where the duct inlet radius is the characteristic length. The divergence angle is 6°. The duct radius 
is then kept constant and a convergent— divergent nozzle with a throat radius of 0.95 is attached. 
The duct exit radius is 0.98 and its total length is 2.9. The divergent part of the duct ensures that 
the formed shock stays in the inlet region. The overall configuration of the duct ensures that the 
supersonic inflow becomes supersonic at the exit. As the computations will show, a small portion 
of the duct exit flow near its CL becomes subsonic at certain times for the specified inflow 
conditions. This configured duct has also been used by Delery et al. [15] for the Euler equation 
computations of supersonic vortex breakdown to computationally model their experimental setup. 

The NS solver uses a grid of 221 x 51 grid points on two axial planes, where 221 points are in 
the axial direction and 51 points are in the radial direction. In the inlet region up to the 0.74 axial 
station, 100 grid points are used; the other 121 points are used in the remaining part of the duct. 
The grid is also clustered at the CL and the wall. The minimum radial grid size at the CL is 0.002. 
The two axial planes are spaced circumferentially at a prescribed angle so that the aspect ratio of 
the minimum grid size will be < 2. The present grid size and distribution resulted from initial studies 
of their effect on the accuracy of the solution. The upstream Mach number is supersonic and is 
kept at 1.75. The initial profile for the tangential velocity is given by 


v 





( 1 ) 


where U x = 1.74, r m = 0.2 and k c = 0.1. The maximum v/U x (swirl ratio /?) is at r = 0.224. The 
radial velocity w at the initial station is equal to zero, and the radial momentum equation is 
integrated to obtain the initial pressure profile. Finally, the density p is obtained from the definition 
of the speed of sound for the inlet flow. With this compatible set of profiles, the computations are 
carried out accurately in time with At = 0.0025. The wall boundary conditions follow the typical 



Fig. 1. Grid of the configured duct (221 x 51). 


CAF 22 - 4 / 5 — N 


610 


C. H. Liu ei ai 



Fig. 2. Streamlines and Mach contours for swirling flow with transient single-bubble breakdown; with 

M x = 1.75, p = 0.32 and Re = 4,000. 

NS solid boundary conditions. These computations have been carried out on the CRAY 2 at NASA 
Langley Research Center. The CPU time is 28 /is per grid point per iteration for the NS calculations. 

Next, we present the results of the computational study on the effects of Re, exit boundary 
conditions and swirl ratio. 

Effect of Re 

The effect of Re on the vortex-breakdown modes is studied by varying Re between 2000 and 
100,000. The Re is based on the radius of the duct inlet. The swirl ratio is kept fixed at 0.32 and 
the downstream-exit conditions are obtained by extrapolating all of the flow variables from the 
cell centers at the exit. 

For Re = 2000, a shock is captured at the duct inlet, but no vortex breakdown is detected. The 
flow at the exit boundary is supersonic. 
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For Re = 4000, Fig. 2 shows snapshots of the streamlines and Mach contours of the solution. 
For this Re, a single breakdown bubble is seen at t = 5 and is convected downstream as time passes. 
This breakdown bubble is formed during the downstream motion of the inlet shock, which reaches 
its maximum downstream displacement at / = 5. Later on, the shock moves upstream (as is seen 
at t = 8) and the breakdown bubble is convected in the downstream direction. Thereafter, the shock 
stays stationary at the inlet and no vortex-breakdown bubbles are formed behind the shock. This 
swirling flow case shows a transient single-bubble breakdown flow. 

For Re = 20,000, Fig. 3 shows snapshots of the streamlines and Mach contours of the solution. 
These snapshots show a vortex-breakdown mechanism of evolution, convection, merging and 
shedding of bubbles, and the inlet shock first moves downstream, then upstream and finally 
downstream. Thereafter, the inlet shock becomes stationary, and no bubbles are formed behind 
the shock. This swirling flow case shows a transient multibubble breakdown flow. 

For Re = 100,000, Fig. 4 shows snapshots of the streamlines and Mach contours of the solution. 
The streamline snapshots show multibubble vortex-breakdown evolution, convection, merging and 
shedding. The time-accurate integration was carried out up to t = 200, and the solution showed 
periodic multifrequency cycles of vortex-breakdown bubbles. An example of the merging of 
vortex-breakdown bubbles of the same sign of vorticity is shown at t = 17. An example of the 



Fig. 5. Enlargement of streamlines of periodic multifrequency, multibubble breakdown; with t = 84 

and 87. 
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convection and shedding of vortex-breakdown bubbles is shown at t = 25. In a comparison of the 
streamline solutions at / = 25 and 89, the solutions are almost identical. This result conclusively 
shows that the breakdown process is periodic. The Mach contours show the dynamics of inlet shock 
motion. In the time range 3 / < 8, the inlet shock moves upstream toward the inlet, and its central 
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Fig. 6. Streamlines and Mach contours for swirling flow with periodic multibubble, multifrequency 
vortex breakdown; with p b = 2p x , Riemann-invariant boundary conditions, M x = 1.75, =0.32 and 

Re = 100,000. 
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portion moves outside the inlet section at t = 8. In the time range 8 < / ^ 25, the inlet shock moves 
downstream with corresponding evolution, convection, merging and shedding of breakdown 
bubbles. In the time range 25 < t < 45, the inlet shock maintains its motion in the downstream 
direction at a slower rate than before, while another shock, which is downstream of the inlet shock, 
appears and also moves in the downstream direction. The evolution, convection and shedding 
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slowly continue until t = 66. In the time range 66 < t 78, the downstream shock disappears, and 
a large vortex-breakdown bubble appears and moves upstream. This motion of the bubble is 
accompanied by upstream motion of the inlet shock (/ = 78). Later, the inlet shock again moves 
in the downstream direction, and the process is repeated. An animated movie has been produced 
to show the breakdown modes until a total dimensionless time of t = 200. 
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Figure 5 shows an enlargement of two snapshots of the streamline solutions at t = 84 and 87. 
At t = 84, five vortex-breakdown bubbles spatially alternate in their sign of vortex strength. Six 
stagnation points exist at the axis. At t = 87, seven vortex-breakdown bubbles and seven stagnation 
points can be seen. The figure shows the merging of two bubbles of same sign of vorticity. This 
swirling flow case shows a periodic multifrequency, multibubble breakdown flow. 

Effect of the type of boundary conditions at the exit 

In the present study, the downstream-exit boundary conditions at the subsonic points are replaced 
by using the Riemann-invariant boundary conditions. The Riemann-invariant boundary conditions 
at the subsonic points are applied with three different values of pressure, called the back pressure 
p b : for the first value, p b = p^\ for the second value, p b = 2p„; and for the third value, p b is obtained 
from dpfdx = const. The other four flow variables are extrapolated from the interior cell centers 
at the duct exit. The Re and p are fixed at 100,000 and 0.32, respectively. 

For p b = Pa up to / = 35, the solution is the same as that where the exit boundary conditions 
are extrapolated from the cell centers. (See Fig. 4.) Thereafter, for t > 35, the inlet shock 
continuously moves in the downstream direction and the vortex-breakdown bubbles move ahead 
of the shock. The shock and vortex bubbles are shed and disappear from the duct at advanced 
levels of time. The breakdown mode is termed as a transient multibubble vortex breakdown. The 
shock/vortex-breakdown/bubble system disappears because the back pressure is too low to support 
the inlet shock sufficiently to keep it in the inlet region. Moreover, the Riemann-invariant boundary 
conditions at subsonic points allow the downstream effects to propagate upstream as time increases. 

For p b = 2 p x . Fig. 6 shows snapshots of the streamlines and Mach contours of the solution. A 
comparison of the present solution with the solution in Fig. 4 shows that the two solutions are 
similar, except that the present solution lags that of Fig. 4 in time. The reason for this lag is that the 
back pressure p b of the present case is larger than that shown in Fig. 4. Moreover, the 
Riemann-invariant conditions at subsonic points allow the downstream effects to propagate 
upstream as time increases. The large back pressure, which is felt upstream, supports the inlet shock 
and keeps it in the inlet region. 

For dpfdx = const, Fig. 7 shows snapshots of the streamlines and Mach contours of the solution. 
A comparison of the present solution with the solution in Fig. 4 shows that the two solutions are 
similar until t = 22. Thereafter, for t > 22, the inlet shock continuously moves in the downstream 
direction with the vortex-breakdown bubbles moving ahead of the shock. Again, as in the case 
of p b = p x , the shock and vortex bubbles are shed and disappear from the duct at advanced 
levels of time. The breakdown is a transient multibubble vortex breakdown. The shock/vortex-break- 
down/bubble system disappears because the back pressure obtained from the dpfdx = const condition 
is too low to support the inlet shock and keep it in the inlet region. Moreover, the Riemann-invariant 
conditions at subsonic points allow the downstream effects to propagate upstream as time increases. 

Effect of the swirl ratio ft 

In this section, the effect of the swirl ratio on the vortex-breakdown modes is studied. The 
downstream-exit boundary conditions are obtained by extrapolating the flow-field variables from 
the interior cell centers at the boundary. The Re for all the cases considered is 100,000; ft is varied 
between 0.2 and 0.38. 

3.5 



Fig. 9. Grid of nozzle flow (221 x 51). 
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Fig. 10. Quasi-axisymmetric flow profiles at * = 0 for a supersonic swirling jet from the nozzle, with 

Mj = 3.0, M x = 2.0 and Re = 296,000. 
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For the fi range 0.2-0.27, a shock is captured at the duct inlet and moves slowly in the 
downstream direction in the subsequent time steps until it stops just before the end of the straight 
inlet portion of the duct. Thereafter, it becomes stationary. No vortex breakdown is detected, and 
the flow at the exit boundary remains supersonic. This result shows that as /? is decreased by 15.6% 
from its original value of 0.32 (Fig. 4), vortex breakdown does not develop. 
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Fig. 1 1 continued opposite. 
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Fig. 11. Streamlines and Mach contours for a supersonic swirling jet from the nozzle; with low-frequency, 
almost single-bubble vortex breakdown and extrapolation from interior boundary conditions. 


For the p range 0.28-0.3, a transient single-bubble vortex breakdown is captured that is similar 
to the case in Fig. 2, where Re = 4000 and p = 0.32. This result shows the dominant effect of the 
swirl ratio on the vortex-breakdown development. With a decrease in p of 6.3% from its original 
value of 0.32, the vortex-breakdown mode changes from a periodic, multifrequency multibubble 
breakdown to a transient, single-bubble vortex breakdown. 

For P = 0.38, Fig. 8 shows snapshots of the streamlines and Mach contours of the solution. The 
vortex-breakdown bubbles are larger than those shown in Fig. 4. In the initial time range up to 
/ = 17, the process of evolution, convection, merging and shedding of vortex breakdown bubbles 
continues, and the inlet shock moves first downstream, then upstream and finally downstream. For 
t > 17, a large vortex-breakdown bubble is formed behind the inlet shock, which oscillates with 
very small amplitude around a mean position. The process of evolution, convection, merging and 
shedding of additional small vortex-breakdown bubbles continues, and the large vortex-breakdown 
bubble oscillates with a very small amplitude around a mean position. This vortex-breakdown case 
introduces a new mechanism that is different from those encountered earlier. 

Vortex Breakdown of a Supersonic Flow from a Nozzle 

In this case, a supersonic swirling jet at Mj = 3, which is issued from a nozzle into a supersonic 
uniform flow of = 2, is considered. A grid of 221 x 51 x 2 in the axial, radial and tangential 
directions, respectively, is used. The computational domain in an axial plane has the dimensions 
7 x 3.5 in the axial and radial directions, respectively, where the nozzle exit radius r = 1 . The 
free-stream Re = 296,000. Figure 9 shows the computational domain and a typical grid for this 
external flow case. The grid is clustered at the nozzle exit and at the CL. 

Figure 10 shows the inflow profiles of the axial velocity, swirl velocity, radial velocity, pressure 
and density, which are taken from the experimental data of Ref. [16]. The initial profiles are used 
as quasi-axisymmetric profiles for the present computations. On the cylindrical boundary (side 
wall) of the flow at r = 3.5, free-stream conditions are imposed that correspond to M x = 2. The 
initial conditions in the computational domain are those that correspond to the free-stream 
conditions at M x = 2. The problem is solved with two types of exit boundary conditions at x = 7. 
First an extrapolation of all five variables from the interior cell center is used; then the 
Riemann-invariant boundary conditions are used. 

Extrapolation from interior cell centers 

Figure 1 1 shows snapshots of streamlines and Mach contours of the solution. The streamlines 
show multibubble breakdown at the early levels. These bubbles develop because of the shock system 
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Fig. 12. Streamlines and Mach contours for a supersonic swirling jet from the nozzle with low-frequency, 
almost single-bubble vortex breakdown and Riemann-invariant boundary conditions. 


that is formed at the nozzle exit in the vicinity of the CL. A strong portion of the shock exists at 
the CL, which splits into two oblique shocks; one is a weak shock and the other is a strong shock. 
Vortex-breakdown bubbles develop behind the strong shock. Thereafter, for t > 5, the shock 
system and the vortex-breakdown bubbles move slowly in the downstream direction. At / > 55, 
both the shock system and the vortex-breakdown bubble move upstream. The slow motion of the 
shock system and the vortex-breakdown bubble continues back and forth between these two 
locations. No vortex shedding has been captured during the computations for this case. Most of 
the exit points are continuously supersonic; hence, no upstream effects exist except within a very 
thin layer around the CL. 

Riemann-invariant boundary conditions 

Next, the boundary conditions at the exit are replaced by the Riemann-invariant boundary 
conditions with p b — p x at the subsonic points. Figure 12 shows snapshots of the streamlines and 
Mach contours of the solution. A comparison of the present solution with the previous case shown 
in Fig. 1 1 shows that the present boundary condition has only a slight effect on the solution. This 
result is not unusual because the subsonic region at the exit is very small and, moreover, the exit 
boundary is far from the nozzle exit. 

Figure 13 shows an enlargement of the Mach contours at t = 55 for the flow case shown in 
Fig. 11. The shock system near the nozzle exit is clearly seen. 

CONCLUDING REMARKS 

The numerical simulation and the study of supersonic vortex-breakdown phenomena have been 
examined for internal and external supersonic swirling flows. A time-accurate solution of the 
unsteady, compressible, full NS equations is used to produce the solutions. The equations are 
solved for laminar flows with an implicit, upwind, flux-difference splitting, finite-volume scheme. 
The solutions are obtained for quasi-axisymmetric flows with a three-dimensional code, FTNS3D, 
by forcing the flow-field vector to be equal on two axial planes in close proximity to each other. 
Quasi-axisymmetric flow solutions require 1 order of magnitude less in computational time than 
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Fig. 13. Enlargement of Mach contours at t = 55 for a supersonic swirling jet from the nozzle. 


the three-dimensional flow solutions and still provide substantial physical understanding of the 
supersonic vortex breakdown and the dominant parameters that affect it. 

In the present study, two supersonic swirling flow cases are considered: the first case is a supersonic 
swirling flow in a configured duct; the second case is a supersonic swirling jet flow that is issued 
from a nozzle into another supersonic uniform flow of lower Mach number than the nozzle flow. For 
the configured duct, the effects of Re, the type of downstream-exit boundary conditions and the swirl 
ratio p are studied. As Re is varied from 4000 to 100,000, different modes of vortex breakdown are 
obtained: a transient single-bubble breakdown; a transient multibubble breakdown; and a periodic 
multifrequency, multibubble breakdown. These solutions have been obtained with extrapolated flow 
conditions from the interior cell centers at the exit. For the flow case with Re = 100,000, the 
downstream-exit boundary conditions have been replaced with the Riemann-invariant boundary 
conditions with p b — /? x , p b = 2/7 x and dpjdx = const. The solutions have shown substantially differ- 
ent vortex-breakdown modes which are dependent upon the type of exit boundary conditions. The 
reason for this result is the upstream effect of the type of exit boundary condition at the exit subsonic 
points. Again, for the flow case with Re = 100,000, p has been varied from 0.2 to 0.38. No vortex 
breakdown develops in the p range 0.2-0.27. In the P range 0.28-0.3, a transient single-bubble 
breakdown develops. At p = 0.38, a quasi-steady, large vortex-breakdown bubble develops with 
small bubbles that experience convection, merging and shedding around the large bubble. 

For nozzle jet flow, the type of downstream-exit-boundary condition has very little effect on 
the vortex-breakdown mode. This result occurs for two reasons: first, most of the exit portion of 
the flow is supersonic, and only a very thin subsonic portion exists around the CL, second, the 
downstream exit is located at a large distance from the nozzle exit therefore, the upstream 
propagation of the type of exit boundary condition at the subsonic points is very small. 
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Abstract. Steady and unsteady asymmetric vortical flows around slender bodies at high angles of attack are 
solved using the unsteady, compressible, this-layer Navier-Stokes equations. An implicit, upwind-biased, 
flux-difference splitting, finite-volume scheme is used for the numerical computations. For supersonic flows 
past point cones, the locally conical flow assumption has been used for efficient computational studies of this 
phenomenon. Asymmetric flows past a 5° semiapex-angle circular cone at different angles of attack, free-stream 
Mach numbers, and Reynolds numbers has been studied in responses to different sources of disturbances. The 
effects of grid fineness and computational domain size have also been investigated. Next, the responses of 
three-dimensional supersonic asymmetric flow around a 5° circular cone at different angles of attack and 
Reynolds numbers to short-duration sideslip disturbances are presented. The results show that flow asymmetry 
becomes stronger as the Reynolds number and angles of attack are increased. The asymmetric solutions show 
spatial vortex shedding which is qualitatively similar to the temporal vortex shedding of the unsteady locally 
conical flow. A cylindrical afterbody is also added to the same cone to study the effect of a cylindrical part on 
the flow asymmetry. One of the cases of flow over a cone-cylinder configuration is validated fairly well by 
experimental data. 


1. Introduction 

Most flight vehicles are designed for attached flow at low angle-of-attack cruise conditions. 
However, for fighter aircraft or missiles under maneuvering conditions, the high angle-of-at- 
tack flight regime is of vital importance. At high angle of attack, slender bodies and highly 
swept wings, common to both fighter aircraft and missiles, led to extensive regions of vortical 
flow on the leeside of the body because of three-dimensional boundary-layer separation. If 
the vortices are both symmetric and stable, their influences can be exploited favorably to 
provide high lift and maneuverability for the vehicle. The region of favorable influence is 
terminated by the onset of asymmetric vortices and the occurrence of vortex breakdown. Such 
phenomena produce large side forces and moments, which may be larger than those 
attainable by the vehicle control system, thus jeopardizing flight safety. 

In the next section, the physical characteristics of vortical flows about various slender 
bodies are described. This is followed by a survey of the experimental and computational 
research work on asymmetric vortex flows. 
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1.1. Physics and experimental work 

Keener and Chapman (1977) categorized vortical flow regimes into four distinct flow 
patterns for slender bodies at various angles of attack (with zero sideslip). These patterns also 
reflect the diminishing influence of the axial flow component (fig. 1). The first pattern 
develops in the very low angle-of-attack range, where the flow is attached and vortex free, and 
the axial flow is dominant. At moderate to high angles of attack, the crossflow influence 
becomes of the same order of magnitude as that of the axial flow, and large scale vortices are 
formed on the leeward side of bodies because of three-dimensional boundary-layer separa- 
tion. In this angle-of-attack range, the vortices are both stable and symmetric, and the large 
increments in normal force due to the low pressure induced on the leeward surface by the 
vortices can be exploited to aerodynamic advantage. A majority of the research work in 
vortical flows has been focused on understanding this symmetric flow pattern. At even higher 
angles of attack, the crossflow effects start to dominate and the vortices may lose their 
stability or even symmetry, which may lead to asymmetric vortices about a symmetric body or 
breakdown of the vortices. Either phenomenon may occur in a quasi-steady or unsteady 
fashion. Both the asymmetric disposition of the vortices and vortex breakdown give rise to 
sudden and potentially catastrophic changes in side-force and moment characteristics. Hence, 
prediction and understanding of the onset of vortex asymmetry and vortex breakdown are 
essential. The fourth flow pattern develops at extremely high angles of attack (up to 90°), 
where the crossflow influence dominates completely, and the leeside flow is characterized by 
an unsteady diffuse wake, with the possibility of having either random or periodic vortex 
shedding depending upon the Reynolds number, Mach number, and geometric details. The 
asymmetric time-dependent vortex shedding is similar to the von Karman vortex sheet in 
two-dimensional flows around cylinders. 

Historically, highly swept, round and sharp leading-edge wings and pointed slender bodies 
are common generic models for the principal components of real fighter aircraft and missiles. 
The study of vortical flows around these isolated aerodynamic components plays an important 
role in the understanding of vortex flows under various conditions including unsteady 
vortex-dominated flows, vortex/shock interaction, asymmetric vortex flow, and vortex break- 
down. For the design of modern fighter aircraft and missiles, the prediction of the onset of 
vortical flow asymmetry is essential. For isolated pointed forebodies, the onset of asymmetry 
occurs when the relative incidence (ratio of angle of attack to semiapex angle of the forebody) 
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Fig. 1. Effect of angle of attack on leeside flow field. 
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exceeds a certain value, e.g., for a pointed circular cone, the relative incidence must be higher 
than two, which has been documented by Peake and Tobak (1982). This flow pattern develops 
about symmetric slender bodies at zero-degree sideslip in response to small perturbations in 
body geometry at the nose or in the flight conditions such as transient sideslip and acoustic 
disturbances. The sudden changes in side force and moment characteristics resulting from the 
asymmetry, in many instances, are sufficiently large to trigger an aircraft or missile to spin. At 
relative incidences near the onset of the asymmetry, the flow is nominally steady. At 
sufficiently high relative incidences, the flow becomes unsteady and asymmetric, with vortex 
shedding either randomly or periodically. 

The literature and recent research work, both computational and experimental, show 
extensive work in the area of study of symmetric vortex flows. Surprisingly, very limited 
research work exists in the area of steady and unsteady asymmetric flows. Recently, a small 
number of computational research studies by several investigators have focused on predicting 
and analyzing the onset of flow asymmetry over slender bodies. This asymmetric vortex 
formation is still an outstanding problem whose physics are poorly understood. While 
experimental studies have produced flow visualization of steady and unsteady asymmetric 
flows on slender bodies, the mechanisms which lead to flow asymmetry are not well 
understood. 

Currently, two mechanisms exist in the literature for explaining the evolution of asymmetry 
(for example. Peak and Tobak, 1982; Skow and Peake, 1982; Lamont, 1982; Yanta and 
Wardlaw, 1982). The first of two these hypotheses appears to operate in both the laminar and 
fully turbulent separation regimes. It suggests that the asymmetry occurs because of the 
instability of the velocity profiles in the vicinity of the saddle point that exists in the crossflow 
planes above the projections of the body vortices. The second hypothesis relates the asymme- 
try to the occurrence of asymmetric boundary-layer transition, leading to an effectively 
asymmetric mean flow about a given body. The onset of asymmetry over slender bodies is 
accompanied by a rapid, local asymmetric movement of the secondary separation line and 
then the primary separation lines circumferentially, precipitated by an asymmetric transition 
region. Although the second mechanism is operable only within the transition zone, the 
former mechanism plays a role in both laminar and fully turbulent flows. For pointed slender 
bodies, the first mechanism produces higher side forces than those produced by the second 
mechanism. Indeed, the implications from the experimental work of Lamont (1980, 1982) with 
tangent-ogive cylinders is that the vortex wake is less structured in the transition domain, 
leading to reduced side and normal forces. In the laminar or fully turbulent regions, the 
vortex structure is well organized, giving rise to larger forces. 

The asymmetric vortex wake usually develops from asymmetric separation line positions on 
the body, but the latter does not appear to be a necessary condition for the former to occur. 
Asymmetric flow has been documented for sharp-edge delta wings where the primary 
separation is fixed at the leading edge (for example. Shanks, 1963; Keener and Chapman, 
1977; Ayoub, 1987; Rediniotis and Telionis, 1989). Generally, even though the separation 
lines are fixed at the sharp leading edges, asymmetry occurs at higher relative incidences than 
those obtained with smooth pointed forebodies or forebody-cylinder configurations. The 
occurrence of asymmetry is attributed to the hydrodynamic instability in the vortex flowfield 
resulting from the crowding together of the vortices as the wing semi-nose angle is decreased. 

The obvious challenges to computational fluid dynamicists is to simulate the asymmetric 
vortex flows through the existing two hypothesized mechanisms, which has been discussed 
earlier. The second challenge is to investigate the determinable parameters for the onset of 
vortical flow asymmetry. These challenges represent the motivation behind the present paper. 
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1.2. Computational work 

Recently, several attempts have been carried out to computationally simulate steady and 
unsteady asymmetric vortical flows around slender bodies of revolution. Graham and Hankey 
(1982) presented the first attempt to compute asymmetric vortical flow over a cone-cylinder 
body, which had been tested experimentally by Thomson and Morrison (1971). They used the 
MacCormack explicit finite-difference scheme to solve the unsteady full Navier-Stokes 
equations for a laminar flow on a relatively coarse grid. The computed asymmetric vortical 
flow was found to be numerically induced by the MacCormack algorithm by its noncentered 
spatial differencing. It was believed that a very small perturbation was induced by the 
finite-difference algorithm truncation error, which triggered an instability at the saddle point 
above the body. Hence, the instability was induced by numerical bias which was physically 
amplified to produce the asymmetry. By switching the algorithm’s sweep direction, the 
asymmetry pattern was reversed. Discrepancies between numerical and wind-tunnel results 
were attributed to insufficient grid resolution, since small disturbances would not be amplified 
on the coarse grid. 

In an attempt to simulate asymmetric vortex flow around an ogive-cylinder body at very 
high angles of attack and at subsonic speeds, Degani and Schiff (1989) obtained asymmetric 
flow solutions to the thin-layer Navier-Stokes equations by introducing a forced asymmetric 
disturbance near the body nose in the form of a small surface jet. When the jet was turned 
off, the flow asymmetry was dissipated and the flow recovered its symmetry. In a later paper 
by Schiff, Degani and Gavali (1989), the unsteady, thin-layer Navier-Stokes equations were 
used to compute the same problem. Vortex unsteadiness developed with increasing angles of 
attack. The behavior of the fluctuations with incidence paralleled the trends observed in 
experiments by Degani and Zilliac (1988). Degani (1990) used the same computational 
scheme to predict the flow around the same ogive-cylinder body for angles of attack a = 20° 
to 80°. His numerical experiments were focused on investigating the origin of the vortex 
asymmetry. Based on his results, the flow field around slender bodies was divided into three 
main groups, depending on the angle-of-attack range. In the range 0° < a < 30°, the results 
show that the flow was symmetric and introduction of small disturbances near the nose had 
only a small effect on the flow asymmetry. In the second range, 30° < a < 60°, the flow 
became steady asymmetric upon introduction of a spaced-fixed forced disturbance near the 
nose. However, when the disturbance was removed, the flow recovered its symmetric shape. 
The origin of asymmetry was attributed to a convective-type instability mechanism. In the very 
high range, 60° < a < 80°, the flow became unsteady with vortex shedding upon introduction 
of a small transient disturbance with short duration. The origin of flow unsteadiness and 
vortex shedding was attributed to an absolute-type instability mechanism. Although this 
investigation revealed good tentative conclusions, there are several remaining questions to be 
addressed. These questions are related to the dissipative effects of the scheme, particularly in 
the crossflow planes, and to the grid fineness and its resolution of the disturbance growth. 

Steady solutions of the incompressible, full Navier-Stokes equations for vortical flow over 
a sideslipping delta wing have first been presented by Hsu and Liu (1990). Results were 
compared with measured data for force and moment coefficients as well as vortex-core 
positions. However, the vortical strength was underpredicted, because of either a lack of grid 
resolution in the vortical region or an inadequate turbulence model for this massively 
separated flow. Strong flow asymmetry was obtained due to the 12° sideslip angle. 

Asymmetric vortical flow simulation due to various types of short-duration disturbances 
was attempted by several investigators. Siclari and Marconi (1989) also used the unsteady, full 
Navier-Stokes equations with a multi-grid, central-difference, finite-volume scheme to solve 
for steady, asymmetric, iocally-conical flows around a 5° semiapex-angle cone over a wide 
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range of angles of attack. At very high incidences, a steady asymmetric solution was captured. 
Later, the same scheme was applied to solve for steady, asymmetric, locally conical flows 
around cones with elliptic, diamond, and biparabolic sections (Siclari, 1990). 

The present authors (Liu et al., 1990; Kandil et ai., 1990a-c) investigated the prediction 
and control of asymmetric supersonic vortex flows around circular and noncircular cones over 
wide ranges of angles of attack, Mach numbers, and Reynolds numbers with locally conical 
flow assumptions. Unsteady asymmetric vortex flows with periodic vortex shedding were 
captured using several different schemes. Later, Kandil et al. (1991a) compared the asymmet- 
ric flow solutions using thin-layer Navier-Stokes and full Navier-Stokes equations. The 
three-dimensional asymmetric flow solutions around circular cones and cone-cylinder configu- 
rations were also studied by the Kandil et al. (1991b). A comprehensive review of the research 
work done by the present authors is presented in the next section. 

1.3. Present work 

In this work, the unsteady, compressible, thin-layer Navier-Stokes equations are used to 
study supersonic, asymmetric, vortical flows. The onset of flow asymmetry occurs when the 
relative incidence of pointed forebodies exceeds certain critical values. At these critical values 
of relative incidence, flow asymmetry develops due to natural and/or forced disturbances. In 
actual flows, the origin of natural disturbances may be a transient sideslip, an acoustic 
disturbance, or similar disturbances of short duration. The origin of forced disturbances may 
be geometric imperfections in the nose region or similar disturbances of a permanent nature. 
The present work is focused on the evolution of flow asymmetry due to assumed natural-type 
disturbances. Two types of flow disturbances are studied: a random round-off error disturb- 
ance and a controlled transient sideslip disturbance with short duration. In addition to 
relative incidence as one of the determinable parameters for the onset of flow asymmetry, the 
effects of free-stream Mach number, Reynolds number, and cylindrical afterbody are studied 
and have been determined to be important parameters. 

Because of the expensive computational resources required for solving three dimensional 
problems, the first part of the computational studies have been applied to supersonic, locally 
conical flows around point cones. Therefore, the mechanism for the onset of steady and 
unsteady flow asymmetry can be studied efficiently and delineated by solving the locally 
conical problems before the three-dimensional problems are examined. In the second part, 
three-dimensional asymmetric supersonic flows over a cone and cone-cylinder configurations 
are investigated, based on the study of the locally conical flow solutions. 


2. Formulation 

In high Reynolds number viscous flows the effects of viscosity are mostly concentrated in 
narrow regions adjacent to solid bodies and in narrow regions of freeshear layers. Owing to 
computer memory limitations, only a limited number of grid points is available for clustering 
mesh points in these regions. As a result, fine-grid spacing is used in directions which are 
nearly normal to these regions, and coarse-grid spacing must be used tangent to these regions. 
In boundary-layer theory, perturbation analysis shows that streamwise components of the 
viscous terms can be neglected relative to those in the normal direction. Similar arguments 
can be applied to the Navier-Stokes equations as a justification for the thin-layer approxima- 
tion. The thin-layer approximation is not the same as the boundary-layer approximation, since 
an approximate form of the normal momentum equation is retained and pressure variation 
across the boundary-layer thickness is taken into consideration. The thin-layer approximation 
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breaks down for low Reynolds numbers and a regions where viscous effects become signifi- 
cant in all directions. Of course, the full Navier-Stokes equations can be incorporated if 
sufficient resolution is provided by the limited grid and if the physical situation warrants it. 
Therefore, in the present work, thin-layer Navier-Stokes equations are chosen to formulate 
two- and three-dimensional flow problems. 

2.1. Thin-layer Navier-Stokes equations 


In many computational applications, the body surface is a boundary of the computational 
domain, and hence the use of body conformal coordinates makes the surface boundary 
condition easy to apply. The transformation of the governing equations from the physical 
Cartesian coordinate system (x v x 2 , x 3 ), to time-independent curvilinear coordinates, 
(£‘> £ 2 , n is given by 


q=J 1 q=J 1 


r = r(x v x 2 ,x 3 ). 

Using the above transformation, the thin-layer Navier-Stokes equations are 
dq/dt + bE m / ar - a( £ v ) 3 /a£ 3 = 0, m = 1,2,3, 
where the flowfield vector, q, is given by 

P 

pu, 

pu 2 
pu 3 

\ e « / 

the inviscid fluxes, E m , are given by 

' pu m ' 
p“iV m + €”p 

pu 2 U m +Z? 2 p 
pu 2 U m + Z”p 

Umi^x+P) j 

and the viscous and heat-conduction flux in £ 3 direction, ( E v ) 3 , is given by 

0 


E m -r 


(^W 1 


tx, T jl 

fyj . 3 

l^J 


; = 1,2,3. 


(1) 

( 2 ) 


(3) 


(4) 


(5) 


The contravariant velocity component in the £ m direction is 
U =£ m u 

Sxj U J> 

and any element corresponding to the three momentum equations in eq. (5) is given by 


( 6 ) 




Re 


+ /=123 
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where f * = 9£ 3 /9;t ; . The last element in eq. (5) is 


= 


Re 


9u 


3 u L 


7 I/3 + £ 3 I T + 

; 9£ 3 3 b * a ^ 3 


1 


3n 2 


^ (y-l)Pr df 


( 8 ) 


In eq. (5), the ts represent the Cartesian components of the shear-stress tensor for a 
Newtonian fluid, assuming Stokes hypothesis and the bs are the shear-dissipation power and 
conduction heat transfer. The inverse of the Jacobian matrix of transformation is 


3(*i, * 2 > * 3 ) 

W7K¥) 


and the metric terms are 



X H 2 

X U 

*2 (' 

X 2 f2 

X 2 ( 

X 3!' 

x ie- 

X K 


t _ bx, bx k 

^- = 2 ^ e Un e lk „— - — , 


(9) 


( 10 ) 


where e ijn and e ikm are the permutation symbols. 

For convenience, all the tensors are expressed in indicial notation. The flow variables are 
introduced in non-dimensional form, and each is referenced to its appropriate free-stream 
value. The non-dimensional density, p, Cartesian velocity components u u u 2 , u v total energy 
e t , viscosity p, and speed of sound a, are defined as the ratio of the corresponding physical 
quantities to those in the free stream, namely p a x , p x al, p^, and a x , respectively. The 
pressure, p, is non-dimensionalized by p x a*, and is related to the total energy for an ideal gas 
by the equation 


p = (y - l)(e, - jpUjUj), (11) 

where y is the ratio of specific heats, and its value is taken to be 1.4 in the present research 
work. The coordinates jc,, x 2 , x 3 , and time, f, are non-dimensionalized by a characteristic 
length, L, and a characteristic time, L/a x , respectively. The viscosity, p, is evaluated by using 
Sutherland’s law 


p = T 3/2 [(l + C)/(T+ C)], (12) 

where T is the temperature and C is the Sutherland constant, which is 110.4 K. The Prandtl 
number, Pr, is chosen to be 0.72. The Reynolds number is defined as Re = pJJ^L / p x , and 
the characteristic length, L is chosen as the length of the body. 

The values of all the free-stream flow quantities which are used as the initial conditions for 
all applications are given as follows: 

P : c = 1, U^ = M, cos a cos f3, u 2oc = —M x sin /3, 

« 3 cc = sin a cos /3, e lx = l/y(y - 1) + \M*, = 1/y, 

a x = T x =l, U 0 =(«^) ,/2 , M„ = UJa„, (13) 

where M K is the free-stream Mach number, a is the angle of attack, and /3 is the sideslip 
angle. 


22 Locally conical Navier-Stokes equations 

For supersonic flows, the three-dimensional Navier-Stokes equations can be transformed 
into the simpler conical flow equations by using the conical coordinates, X , Y, and Z, with 

X = x x , Y = x 2 /x ,, Z=x 3 /x x . (14) 
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Physically, a conical flow has the property that all flow quantities are invariant along rays that 
emanate from the apex of the conical surface. Using eq. (14) to transform the full Navier- 
Stokes equations to A, Y, Z coordinates and imposing the conical flow property, the resulting 
equations in abstract form are given by 


3 q 3(F-F V ) 3(G — G v ) 

*¥ + »F^ + ' L ir _L + 2(£l ~ £ - ,)_0, 


where the inviscid fluxes are 
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and the viscous fluxes are 
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(19) 


The shear stresses, dissipation power, and heat transfer terms are obtained by using chain-rule 
differentiation and enforcing the conical flow property, i.e. all derivatives in the A-direction 
are zero. For example, the principal stress, r n , can be simplified as 

„^ U 2 bu 2 dt M 

1 2 Y — + 2 Z — H 1- — I . 

\ ay a z ay a z ) 


T n = - 


Re A 


( 20 ) 


The resulting equations (15) have spatial variation in the Y- and Z-directions only. Thus, 
these equations are two-dimensional equations with source terms. Hence, they are more 
economical to solve than the three-dimensional equations. It is also noticed that the time-de- 
rivative term in eq. (15) is multiplied by X and an axial length-scale dependence exists in the 
viscous terms (20). Hence, eq. (15) is not self-similar in A-direction, and thus it does not 
represent a globally conical flow. Only the steady inviscid flow equation represents a globally 
conical flow. However, for unsteady viscous flow over a conical body, if A is fixed at a certain 
location, the flow may be thought of as “locally conical”, with the Reynolds number 
determining the location of the conical plane in which eq. (15) is solved. The best that can be 
done to make use of this equation is to select a constant value for A, and solve the resulting 
equation for what we call “locally conical flow”. 
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3. Computational scheme 


The computational scheme used to solve the governing equations is an implicit, upwind-bi- 
ased flux-difference splitting, finite-volume scheme. The discretized equation is integrated 
numerically in time using the Euler implicit time-differencing method. The linearized, 
backward-time approximation for the flowfield vector is written in the delta form as 


where 



Aq n 


R(q n ), 


( 21 ) 


Aq n = q n + 1 - q n , 

«(«")- - [s { , (£,)"- < 22 > 


In eq. (22), 8 ( , is the spatial difference operator. The convective and pressure terms are 
discretized using the flux-difference splitting scheme of Roe and differenced using the 
MUSCL (Monotone Upstream-centered Scheme for Conservation Laws) of van Leer. The 
flux-difference splitting of Roe is based on the approximate Riemann problem. The Riemann 
problem is used as a mechanism to divide the flux difference between the neighboring states, 
such as the interface of two computational cells, into component parts associated with each 
wave field. As each eigenvalue is also associated with its own wave field, so the splitting can 
be done based on the eigenvalues. The smooth flux limiters are used to eliminate oscillations 
in the shock region, and the viscous and heat-flux terms are centrally differenced. The 
resulting difference equations are solved using a spatially split approximate factorization 
along the £ 2 , directions, respectively. The scheme is first-order accurate in time and 
third order accurate in space. Details of the above described scheme are given by Wong 
(1991). 

Since the applications in this paper cover some locally conical flow problems, locally 
conical flow solutions can be obtained by solving the problem in three conical planes using a 
three-dimensional solver. This is achieved by setting the conserved components of the 
flowfield vector, q, to be equal at two planes. All of the locally conical solutions in the present 
work are obtained in this way. 


3.1. Initial and boundary’ conditions 

All the numerical calculations of the steady-flow problems are obtained by using impul- 
sively-started initial conditions, i.e. bodies are suddenly placed in the free stream at angles of 
attack specified by the problem. For unsteady-flow problems, solutions obtained from the 
pseudo time-stepping calculation corresponding to the same flow conditions are used as initial 
conditions in order to save the computational cost for the transient state. 

The boundary conditions for the present work are implemented explicitly. On the solid 
boundary, the no-slip and no-penetration conditions are enforced, i.e. u, = u 2 = « 3 = 0, and 
the normal pressure gradient is assumed to be zero. The adiabatic condition is maintained on 
solid surface. 

To obtain a locally conical flow or three-dimensional solution for supersonic free-stream 
Mach numbers, the computational domain is extended far enough to permit capture of the 
bow-shock formed outside of the body as part of the solution. Since the disturbance from the 
body will not propagate beyond the bow-shock in the crossflow plane, the conditions outside 
the conical shock are the same as the free-stream conditions. Therefore, the farfield boundary 
conditions are specified to be the free-stream conditions. Since the locally conical flow 
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solutions are obtained by solving the problem in three conical planes using a three-dimen- 
sional solver, free-stream conditions are enforced on the farfield boundary of the first and 
third conical planes. For the three-dimensional flows, first-order extrapolation from the 
interior points at the outflow boundary is used. At the geometric plane of symmetry, periodic 
conditions are used since the problem is solved for the whole computational domain. 


4. Results and discussion 

4. 1. Locally conical flows 

The mechanisms which lead to steady and unsteady asymmetric vortical flows past slender 
wings and bodies at high angles of attack at zero sideslip are not well understood. From 
experimental studies of these phenomena, several investigators proposed two mechanisms for 
explaining the origin of the flow asymmetry. These have already been described in the 
previous section. The first mechanism, asymmetric flow due to a saddle-point instability, is 
demonstrated in this section. Two types of flow disturbance, a random round-off error 
disturbance and a controlled transient sideslip disturbance with short duration, are used to 
demonstrate the mechanism which leads to flow asymmetry. In addition to the relative 
incidence as one of the determinable parameters for the onset of flow asymmetry, other 
influential parameters such as the Mach number are studied and presented in this section. 

4.1.1. Steady asymmetric flows over a cone 

Supersonic flows over a 5° semiapex angle cone at a Reynolds number Re = 10 5 have been 
computed. The grids used in all the numerical tests in this section are generated by using the 
modified Joukowski transformation with a geometric series for grid clustering near the cone 
surface. For all the cases, a grid of 161 X 81 points is used, where the first number is the 
number of points around the cone and the second number is the number of points normal to 
the cone surface. A 241 X 121 grid and a 161 X 81 grid with different mesh fineness ratios or 
different computational domain sizes have also been used to test the effect of grid fineness 
and domain size on the numerical solutions. A typical grid of 161 x 81 points is shown in 
fig. 2. 
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total-pressure-loss contours 




surface-pressure coefficient 

(a) = 10- (b) A£,,„ = 10- (c) = 10_S 

Fig. 3. Effect of minimum spacing A£^ in on the asymmetric solution (a = 20°, A/* =1.8, Re = 10 s , 161x81 

points,/?/ = 21r). 
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To establish an optimum grid and to ensure that the asymmetric flow solution is unique, 
irrespective of the grid fineness and the computational domain size, numerical tests have been 
carried out with several different grids. The tests have computed the supersonic flow around a 
cone at a = 20°, = 1.8, and Re = 10 5 (relative incidence is four for this case). A grid of 

161 x 81 points in the circumferential and normal directions, respectively, has been used with 
different minimum grid spacing, A£ 3 , at the solid boundary, while the maximum radius of the 
computational domain, r h is fixed at 21 r, where r is the radius of the circular cone at the 
axial station of unity. Three cases, computed using A£ 3 = 10 -3 , 10 -4 , and 10 _5 , are shown in 
fig. 3. In fig. 3, the logarithmic residual error versus the number of iterations, the surface 
pressure versus the azimuthal angle, 0, which is measured from the leeward plane of 
geometric symmetry, and the total-pressure-loss contours are shown. The residual error 




convergent history surface-pressure coefficient 



total-pressure-loss contours 

Fig. 4. Effect of increased grid density on the asymmetric solution (a = 20°, = 1 . 8 . Re = 10 5 , 241x121 points. 

A^min “ 10 -6 , R t = 21 r ). 
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figures show that the error reaches machine zero (10 — 10 to 10" n ) in about 2500 time steps in 
all cases and the solutions are symmetric at this point. Afterwards, with the machine 
round-off error is acting as a random disturbance to the flow field, the residual error grows, 
then drops down by at least another seven orders of magnitude, and finally stays constant 
thereafter (constant residual error for 2000 iterations is shown). The pressure coefficient and 
total-pressure-loss contours show that the flow becomes steady, asymmetric, and stable. The 
solution of the three cases are not necessarily the same because the source of disturbance is a 
random one, and it is possible that the solutions are mirror images of each other. Other types 
of disturbances will be discussed in the next section. Furthermore, a grid of 241 x 121 points 
with minimum spacing of 10 -6 is used to test the effect of grid density on the asymmetric 
solution. Figure 4 shows the results of this case. The residual error figure shows that the error 
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Fig. 5. Effect of increased computational domain size on the asymmetric solution (a = 20°, A/ x =1.8, Re = 10\ 

161x81 points, A^ in = 10' 4 , Rj = 32 r). 
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drops ten orders of magnitude in 7500 time steps, then grows about five orders of magnitude 
after being triggered by the machine round-off error and then converges to the same 
asymmetric solution. Since the asymmetric solution is unique, irrespective of the size of 
minimum grid spacing and grid density, an optimum grid spacing of 10 ~ 4 is chosen in the 
present study. 

Two grids of 161 X 81 points with the maximum computational domain radius increased 
from 21 r to 32 r are used to test the effect of the domain size on the solution. The optimum 
minimum spacing is used for the grid. The results of this case are shown in fig. 5. The residual 
history shows a similar trend in going through a symmetric unstable solution and then to an 
asymmetric stable solution. The pressure coefficient and total-pressure-loss contours figures 
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Fig. 6. Effect of spatial disturbance on the asymmetric solution (a = 20°, Af„=1.8, Re = 10 5 , 161x81 points, 

Afita-10- 4 , */-21 r>. 



C.H. Liu et al. / Vortical flows around slender bodies 


423 


are consistent with the results shown in fig. 3b. Thus, the optimum grid spacing of A£ 3 = 10 
and the maximum radius of 21 r are chosen to be used for all the cases presented in this 
paper. 

As mentioned in the previous section, the locally conical flow solution is obtained by 
forcing the equality of the flow-field vector at two cross sections, which are taken as £ 1 = 0.95 
and 1. A numerical test has been performed for the same flow conditions except that the 
solution is achieved by forcing the equality of the vector, q, at = 0.995 and 1. The purpose 
of this task is to test the spatial disturbance on the asymmetric solution. Figure 6 shows the 
results of convergence history, pressure coefficient, and total-pressure-loss contours. The 
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Fig. 7. Steady asymmetric flow solution for circular cone due to ±0.5° transient sideslip (/3) (a = 20°, M»= 1.8, 

Re = 10 5 ). 
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residual error shows that the solution takes larger number of time steps for the asymmetry to 
be triggered by the machine round-off error. However, the surface pressure and the total- 
pressure-loss contours confirm the uniqueness of the asymmetric solution. Since A£' = 0.005 
is a small disturbance to the locally conical flow assumption, it is reasonable to have longer 
time steps to obtain the asymmetric solutions. To efficiently use of the limited computational 
resources, A^ 1 = 0.05 is used for all the locally conical flow problems in the present work. 

Since the magnitude of residual errors shown in the above cases is so small, it is believed 
that the disturbance which triggered the flow asymmetry can be attributed to the machine 
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Fig. 9. Comparison of crossflow velocity vectors and total-pressure-loss contours for circular cone at different Mach 

numbers (a - 20°, Re = 10 5 ). 
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round-off error. This type of disturbance is random in nature. In the next section, the results 
of a controlled type of disturbance are presented. 

4.1.2. Controlled transient sideslip disturbance 

In this section, steady asymmetric flow solutions due to a transient sideslip disturbance of 
short duration are presented. Results of the transient sideslip, fi = ± 0.5°, are shown in fig. 7. 
The residual error figures show a drop of seven orders of magnitude in the first 2000 time 
steps. At this step, a sideslip disturbance is imposed for six time steps, then it is removed. 
Irrespective of the magnitude or the direction of the sideslip disturbance, the residual error 
increases by six orders of magnitude, then drops down very rapidly and a stable asymmetric 
flow solution is obtained. The asymmetric solutions corresponding to the ±0.5° sideslip 
disturbances are mirror images of each other, as can be seen from the surface-pressure 
distributions, crossflow velocity vectors, and total-pressure-loss contours. Moreover, the final 
stable asymmetric solutions of the ± 0.5° sideslip disturbances are the same or mirror images 
as those from random disturbances shown in figs. 3-5. 

4.1.3. Steady asymmetric flow at different Mach numbers 

Using the same optimum grid and the same 5° semiapex angle cone at a = 20°, three cases 
of locally conical flow solutions with free-stream Mach numbers ranging from 2.2 to 3.0 have 
been computed. The effect of the free-stream Mach number on the convergence history, 
surface pressure, crossflow velocity, and total-pressure-loss contours are shown in fig. 8 and 9. 
At M „ = 2.2, the residual error shows that the stable asymmetric flow is obtained within the 
same number of time steps as that of the M m = 1.8 case. At = 2.6, the residual error shows 
that the final asymmetric solution is obtained after a larger number of time steps. At 
M„ = 3.0, no asymmetric flow has been captured and the flow stayed symmetrically stable. 
The surface pressure figures show that the flow asymmetry gets weaker as the Mach number 
is increased. This conclusion is strongly supported by the crossflow velocity vectors and the 
total-pressure-loss contours, as shown in fig. 9. It is also noted that since the nature of 
disturbance is random, flow asymmetry changes sides as the Mach number increases, until it 
disappears. The significant feature of these numerical tests is that the asymmetric/symmetric 
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Fig. 10. Time histories of residual error and lift coefficient for unsteady asymmetric flow around circular cone 

(a = 30°, = 1.8, Re = 10 s ). 
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behavior of the solutions is continuous and the general trends with Mach number are in 
agreement with the experimental observations by Peake and Tobak (1982). 

4.1.4. Unsteady asymmetric flows over a cone 

Keeping the Mach number at 1.8 and Reynolds number at 10 5 , the angle of attack is 
increased to 30° for the flow around the same circular cone (relative incidence is six for this 
case). The histories of the logarithmic residual error and the lift coefficient versus the number 
of iteration up to 15 900 time steps are shown in fig. 10. First, pseudo-time stepping has been 
used up to 8000 iterations and the solution has been monitored every 500 iterations. The 
solution is still symmetric at 3000 iterations. Thereafter, the flow asymmetry has been 
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Fig. 12. Snapshots of total-pressure-loss contours for unsteady asymmetric flow around circular cone (a = 30°, 

A/. -1.8, Re = 10 5 , At = 10 -3 ). 
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Fig. 13. (Continued). 


15,700 


obtained through the random disturbance of the scheme. The asymmetry changes randomly 
from the left side to the right side, which indicates a possibility of unsteady asymmetric vortex 
shedding. Therefore, the computations have been restarted from the 8000th time step using 
time-accurate stepping with a minimum global time step, (A f) mjn = 10~ 3 . The residual-error 
and lift-coefficient histories show that, after switching to the time-accurate stepping, a short 



Fig. 14. Snapshots of total-pressure-loss contours for unsteady asymmetric flow around circular cone within one cycle 
(cylinder axis is a time axis; a = 30°, - 1.8, Re = 10 \ A / = 10 -3 ). 
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transient response is followed by a periodic response. Figures 11, 12 and 13 show snapshots of 
the time history of the solution for the surface-pressure coefficient, total-pressure-loss 
contours, and crossflow velocity vectors. The solutions are shown every 100 time steps starting 
from the time step of 15000 to 15 700. At a time step of 15000, the asymmetric flow is seen 
with a vortex already shed from the right side. As time passes, the shed vortex is convected 
into the flow, and the primary vortex on the left side stretches, while the primary vortex on 
the right gets stronger, as seen from the surface-pressure curves in fig. 11. At a time step of 
15 600, the primary vortex on the left side is about to be shed. At the time step of 15 700, the 
primary vortex on the left side is shed into the flow field. It is also noted that the solution at 
the time step of 15 700 is a mirror image to that of the 15000 time step. Hence, the solution 
from the 15 000 to the 15 700 time steps represents the one half cycle of shedding. The 
periodicity of the shedding motion has been captured. The period of oscillation is 10“ 3 X 1400 
= 1.4, which corresponds to a shedding frequency of 4.488 (Strouhal number). Figure 14 
shows a snapshot of the total-pressure-loss contours over one period on a cylinder, with the 
axis of the cylinder representing the time axis. The present unsteady asymmetric flow solution 
has also been obtained exactly by using the flux-vector splitting scheme with the thin-layer 
Navier-Stokes equations and the flux-difference splitting scheme with the full Navier-Stokes 
equations on a finer grid (Kandil et al., 1991). Hence, the present solution is unique and 
independent of the computational scheme or the approximation level of the Navier-Stokes 
equations. 
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4.2. Three-dimensional flows 

Steady and unsteady solutions using the locally conical flow assumptions for supersonic 
flows around a circular cone have been presented in the previous section. For such problems, 
these assumptions reduce the computational time and memory. However, the three-dimen- 
sional effects cannot be neglected in certain vortical regions, such as for flows with massive 
three-dimensional separation, vortex shedding, vortex breakdown, and others. In these re- 
gions, one has to rely on solutions of three-dimensional Navier-Stokes equations. 

In this section, solutions of three-dimensional asymmetric flows around a 5° semiapex-an- 
gle circular cone and two cone-cylinder configurations are presented. Two issues concerning 
the flow asymmetry around a circular-section cone in response to a short duration transient 
sideslip disturbance are addressed. First, for the same cone section and for the same flow 
conditions and disturbance, does the three-dimensional flow solution produce the same 
solution as that of the locally conical solution presented in the previous section? Second, what 
are the effects of angle of attack, Reynolds number, and cylindrical afterbody on flow 
asymmetry? Finally, the results of the asymmetric flow solution are validated with those of the 
experimental data by Landrum (1977). 

4.2. 1. Steady asymmetric flow over a cone 

An O-H grid of 65 X 161 X 81 points in the streamwise (£ l ), circumferential (£ 2 ), and 
normal (£ 3 ) directions, respectively, has been used. The grid is generated in the crossflow 
planes using a modified Joukowski transformation which is applied locally at the grid length 
stations, with algebraic stretching at the cone surface. The crossflow grid (161 X 81) is of the 
same size as that used for the locally conical solutions. In order to retain the same resolution 
for each conical section, the outer boundary is a conical surface with the maximum radius of 
3 L at the cone base, where the L is the length of the cone. The minimum spacing at the cone 
surface ranges from 10" 5 at the cone base to 10“ 6 at the cone apex. In the circumferential 
direction, the grid is equally distributed for the whole computational domain. A typical grid is 
shown in fig. 15. 

For the same flow conditions, a = 20°, = 1.8, and Re = 10 5 , at which the locally conical 

flow solution is asymmetric, a symmetric flow solution has been obtained using the three-di- 
mensional calculation. The difference is explainable as a Reynolds number effect, since the 
locally conical solution is obtained at a fixed axial station, £* = 1.0. As mentioned in the 
formulation, a length scale in the viscous terms (Reynolds number) for steady viscous flow 
remains after the conical transformation. The resulting equations are not self-similar, and the 
location of the conical plane in the transformed equation determines the Reynolds number. 

A slight asymmetric flow solution has been obtained for the three-dimensional cone flow 
after increasing the angle of attack to 40°, and the free-stream Reynolds number to 4 x 10 b , 
and reducing the free-stream March number to 1.4. The flow is assumed fully laminar in the 
numerical computation. During this computation, it has been observed that the computed 
flow remains symmetric about the geometric plane of symmetry at the leeside of the body. 
The symmetry of the solution is then disturbed by introducing a sideslip angle of 2° to the flow 
field for about 100 time steps and then it is removed. Thereafter, the pseudo-time stepping is 
continued until the residual error drops again four orders of magnitude and a stable 
asymmetric solution is obtained. The total-pressure-loss contours of this case are shown in fig. 
16. Although the resulting flow is no longer symmetric, the asymmetry is relatively small. In 
this case, the vortices still lie close to the leeward-body surface, and the size of the shear layer 
and height of the primary vortices grow with increasing distance downstream. It is also seen 
that the solution is almost self-similar over a long distance of the cone length. 

Next, the Reynolds number is increased to 5 x 10 6 and 6 x 10 6 , keeping the other flow 
conditions the same as those of the previous flow case. Again, the source of the disturbance to 
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Fig. 16. Total-pressure-loss contours of flow around at 5° semiapex-angle cone ( a = 40°, M = 1.4. Re = 4x 10*). 



Fig. 17. Total-pressure-loss contours of flow around 5° semiapex-angle cone (a = 40°, M x - 1.4, Re = 5x 10 6 ). 
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break the symmetry of the solution is a 2° transient sideslip of short duration. The computa- 
tion has been monitored every 300 time steps until there are no significant changes in the 
steady-state solutions. Figures 17 and 18 show the total-pressure-loss solutions for these cases. 
Figure 17 shows that the asymmetry of the vortical flow becomes strong and the self similarity 
of the flow asymmetry is substantially lost. However, it is seen that the primary vortices do not 
change sides as moving in the downstream direction. Figure 18, which corresponds to the 
free-stream Reynolds number of 6 x 10 6 , shows flow asymmetry that changes sides as the 
solutions develop along the length of the cone and as vortices are shed into the flow field. 
Since the solution is steady, the vortex shedding is a spatial one. A close study of the solutions 
shown between the fourth crossflow plane and seventh crossflow plane reveals that the vortex 
shedding changes from the right side (looking downstream, fourth crossflow plane) to the left 
side (seventh crossflow plane). The solutions on these two planes are nearly scaled mirror 
images of each other. The present spatial flow asymmetry is qualitatively similar to the 
temporal flow asymmetry of the locally conical flow solution shown in fig. 14. It is seen that 
most of the asymmetric vortex flow characteristics and physics in the crossflow can be 
captured quantitatively by the locally conical flow solutions. 

Figure 19 shows the total-pressure-loss solution for the same cone for Re = 8 x 10 6 . The 
asymmetry of the vortex flow becomes much stronger, as compared with the previous cases of 
figs. 16-18. By comparing the solution of this case with that of the Re = 6 X 10 6 , it is noticed 
that the flow asymmetry of the case with Re = 8 x 10 6 changes sides along a shorter axial 
distance (third and fifth crossflow planes). Moreover, the flow asymmetry of the case with 
higher Reynolds number changes sides one more time (fifth and ninth crossflow planes) and, 
thus, a complete wave length of flow asymmetry is formed between the third and ninth 
crossflow planes. A close study of the mechanism of spatial vortex shedding along the cone 
reveals that it is similar to the unsteady vortex shedding of the locally conical flow solution. At 
the third crossflow plane (x { /L = 0.2), the asymmetric flow is seen with vortex already shed 
from the right side. Moving downstream, the shed vortex is convected into the flow, and the 
shear layer on the right side stretches, while the primary vortex on the left side gets stronger, 
as seen from the surface-pressure curves in fig. 20. At the fifth crossflow plane (x x /L - 0.4), 
the primary vortex on the right side is about to be shed. At the ninth crossflow plane 
(jr,/L = 0.9), the primary vortex on the right side is almost shed in the flow field, while the 
lower part of shear layer on the same side has stretched and shrunk in thickness. It is also 
seen that at the ninth crossflow plane, the flow is approximately a mirror image of that at the 
third crossflow plane. The behavior of the flow asymmetry over one period in fig. 14 is 
qualitatively similar to that of the flow asymmetry over one wave length in fig. 19. 

4.2.2. Unsteady asymmetric flow over a cone 

In this section, solutions of the unsteady supersonic asymmetric flow around the same 
circular cone at an angle of attack of 50° are presented. The free-stream Reynolds number 
and Mach number of this case are 8 X 10 6 and 1.4, respectively. The present flow case has 
been started from the solution obtained for a = 40°, instead of initializing with free-stream 
conditions everywhere. In addition, this steady asymmetric initial condition can be considered 
as the source of disturbance to the flow field, so the use of transient sideslip disturbance is 
not necessary for this case. 

In the computation of locally conical flow problems it has been shown that, once unsteady, 
asymmetric vortex shedding is initiated, the perturbation can be removed. The vortex 
shedding will continue without the need for any further perturbations since the flow is 
unstable. In order to investigate whether the same phenomena exists for the unsteady 
three-dimensional asymmetric flow, the computation has been first done using pseudo-time 
stepping until the residual error drops three orders of magnitude. The flow asymmetry. 
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Fig. 18. Total-pressure-loss contours of flow around 5° semiapex-angle cone (a = 40°, = 1.4, Re = 6x 10 b ). 
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Fig. 19. Total-pressure-loss contours of flow around 5° semiapex-angle cone ( a = 40°, = 1.4, Re = 8x 10 6 ). 
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changing randomly from the left side to right side, has been captured by the pseudo-time 
stepping calculation, which indicates a possibility of unsteady vortex shedding. The computa- 
tion has been continued using time-accurate calculations with a minimum global time step of 
10 -5 . The unsteady structure of the flow at a = 50° is monitored at different time steps to see 
the mechanism of the unsteady vortex shedding and the unsteady behavior of the vortex 
structure. Due to the fine computational grid spacing at the nose region and the cone surface, 
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Fig. 20. (Continued). 


the allowable computational time step is so small that the calculation becomes prohibitively 
expensive for time-marching more than a small part of a shedding cycle. 

Consequently, only small changes are visible in the instantaneous snapshots of the unsteady 
asymmetric flow solutions presented in this section. The flow at time step 10616 is shown by 
the total-pressure-loss contours (fig. 21), surface-pressure coefficients (fig. 23) and enlarge- 
ment of the total-pressure-loss contours (fig. 24). The same quantities for time step 11 816 are 
shown in figs. 22, 25, 26. Figure 21 shows a strongly asymmetric solution with vortex shedding 
changes little with axial location. The shear-layer thickness for each crossflow station extends 
about one and one half times the local diameter of the leeward plane of symmetry, as 
compared with the case of a = 40°. It is evident that all of the three vortices interact with 
each other in a relatively small distance above of the body surface. The blow-up of 
total-pressure-loss contours at the time step of 11 816 (fig. 26) shows that the flow asymmetry 
changes side at the crossflow station of Jt,/L = 0.1, as compared with the same crossflow 
section in fig. 24. Obviously, the total computed time is too short in terms of physical time to 
draw final conclusions for periodic vortex shedding, but the total-pressure-loss contours 
clearly show that the flow is unsteady asymmetric with a possibility of vortex shedding at each 
axial station. 

All of the numerical results have been obtained using either the Cray-2 supercomputer of 
the NASA Langley Research center or the Cray-YMP supercomputer of the NASA Ames 
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Fig. 21. Total-pressure-loss contours of flow around 5° semiapex-angle cone at time-step of 10616 (a « 50°, M 

Re = 8 X 10 6 ). 



Fig. 22. Total-pressure-loss contours of flow around 5° semiapex-angle cone at time step of 11 816 (a — 50°, 

Re = 8 x 10 6 ). 
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Fig. 24. (Continued). 


Research Center. Each of the steady asymmetric flow cases takes about 100 h of CPU time on 
the Cray-2 and 65 h of CPU time on the Cray-YMP computer. The unsteady asymmetric flow 
case takes over 82 h of CPU time on the Cray-YMP computer for the 11 816 computed time 
steps. 


4.2.3. Asymmetric flow over cone-cylinder configurations 

The effect of a cylindrical afterbody on flow asymmetry is investigated by introducing a 
cylindrical afterbody of unit length to the same circular cone. The flow around the resulting 
cone-cylinder configuration is solved with the flow conditions as for a = 40°, M x = 1.4, and 
Re = 4 x 10 6 , which are the same flow conditions of the isolated unit-length cone shown in 
fig. 16. The source of disturbance is the same 2° transient sideslip. The computed total-pres- 
sure-loss contours are shown in fig. 27. It should be noted that slight flow unsteadiness has 
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been observed during the computations. Comparing the results for fig. 16 with fig. 27, it is 
noticed that the flow asymmetry is stronger for the cone-cylinder configuration than that of 
the isolated conical forebody. It should be noted that subsonic flow region does exist inside 
the conical shock surrounding the cone-cylinder configuration, hence, the downstream cylin- 
drical-afterbody boundary has an upstream effect on the flow. There are two reasons for the 
afterbody to increase the flow asymmetry; the first is the increase of the local angle of attack 
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Fig. 25. Surface-pressure coefficient on 5° semiapex-angle cone at time step of 11816 (a = 50°, M„=1.4, 

Re = 8x 10 6 ). 
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Fig. 25. (Continued). 


of the leeward side of the cylinder, and the second is the surface discontinuity at the 
cone— cylinder juncture. Both of these increase the spatial growth of the flow asymmetry. 

Next, a validation flow case is presented by comparing with available experimental data. 
For this purpose, the cone-cylinder configuration of 0.5 : 0.5 (the ratio of the conical forebody 
to the cylindrical afterbody) is used, which was experimentally tested by Landrum (1977). The 
flow conditions for this case are a — 46.1°, = 1.6, and Re = 6.6 x 10 6 . The Reynolds 

number for this case is based on the total body length. The cone semiapex angle is 9.46°. and 
the numerical computation is assumed to be fully laminar. The problem is solved using a grid 
size of 65 x 161 x 81, which has the same resolution in the crossflow plane as the previous 
cases. Figure 28 presents the computed total-pressure-loss contours, which show a relatively 
weak asymmetry at the nose region and a strong spatially growing asymmetry in the 
downstream direction. Figure 29 shows the surface-pressure coefficient along with the 
experimental data, and figs. 30 and 31 show the total-pressure-loss contours and the total 
Mach-number contours in the crossflow planes at the axial stations of 0.075, 0.225. 0.475, and 
0.775. The computed (solid line) and measured (symbol) surface-pressure coefficients are in 
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Fig. 26. Blow-up of the total-pressure-loss contours on 5 C semiapex-angle cone at time step of 11816 (or -50°. 

M x = 1.4, Re = 8X10'’). 
















C.H. Liu et ai / Vortical flows around slender bodies 


445 



Fig. 27. Total-pressure-loss contours of flow around 5° semiapex-angle cone-cylinder configuration (a = 40°. 

M = 1.4, Re = 4Xl0 6 ). 



Fig. 28. Total-pressure-loss contours of flow around 9.46° semiapex-angle cone-cylinder configuration (a = 46.1°, 

A/. - 1.6, Re = 6.6 X 10 6 ). 
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Fig. 29. Surface-pressure coefficient on 9.46° semiapex-angle cone-cylinder configuration (a = 46.1°, M^= 1.6. 

Re = 6.6 X 10 6 ) (experiment: Landrum, 1977). 


good agreement, with the exception of the region close to the leeward plane of symmetry. By 
studying the blow-up of the total-pressure-loss contours (fig. 30), it is seen that a slight flow 
asymmetry starts at jc,/L = 0.075 and spatially grows in the downstream direction. The total 
Mach-number contours (fig. 31) show that the shocks on the primary vortices are asymmetric 
and change sides as moving downstream, as shown by the results at the axial stations of 0.475 
and 0.775. It is noticed that the flow asymmetry of this case is relatively weaker than that of 
the 5° semiapex-angle cone-cylinder case because the relative incidence of the forebody of the 
former case is lower than the latter. 


5. Conclusions 


The main goal of the present work is to predict asymmetric vortex-dominated flows around 
slender bodies over a wide range of angles of attack, Mach numbers, and Reynolds numbers. 
In this section, a summary of the findings of the numerical investigation is presented. First, 
steady and unsteady solutions of supersonic asymmetric flows around a circular cone have 
been obtained using the thin-layer Navier-Stokes equations along with the locally conical flow 
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Fig. 30. Blow-up of the total-pressure-loss contours on 9.46° semiapex-angle cone-cylinder configuration (a = 46.1°. 

= 1 .6. Re = 6.6x10'’). 


assumption. The results have shown that the onset of flow asymmetry occurs when the relative 
incidence of cones exceeds certain critical values. At these critical values of relative incidence, 
asymmetric flow develops, irrespective of the sources of disturbance. Two types of flow 
disturbances of short duration are used to demonstrate that the asymmetric solution is unique 
and that the mechanism leads to flow asymmetry due to instability of the saddle point, even 
without the presence of any permanent disturbance. It has also been shown that as the Mach 
number increases, vortex flow asymmetry becomes weaker. In the high angle-of-attack regime, 
unsteady asymmetric flow with periodic vortex shedding has been uniquely captured. 

Second, asymmetric supersonic three-dimensional flow problems have been solved using 
the thin-layer Navier-Stokes equations. Steady and unsteady flow solutions for asymmetric 
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Fig. 31. Blow-up of the total-Mach-number contours on 9.46° semiapex-angle cone-cylinder configuration (a = 46.1°, 

m x = 1.6. Re = 6.6 x 10 A ). 


supersonic flow over a circular cone and cone-cylinder configurations have been presented. 
However, because there is a serious lack of steady and unsteady three-dimensional detailed 
experimental measurements for supersonic asymmetric vortex flows, only surface-pressure 
coefficients from one of the cases are validated with the experiment in the present study. 

It is shown that the three-dimensional flow calculation does not produce the same solution 
as that of the corresponding flow case under the locally conical assumption. The reason is that 
for the viscous flow problem the transformed equation using the locally conical flow assump- 
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tion is not self similar. However, most of the characteristics and physics of the asymmetric 
flow can be simulated. The three-dimensional numerical results also show that the onset of 
steady and unsteady asymmetric flow develops due to a transient sideslip disturbance of short 
duration, provided that the body is at a certain critical range of Mach number, Reynolds 
number, and angle of attack. The steady and unsteady asymmetric flow solutions have been 
obtained without the need to impose any permanent type disturbance. As the free-stream 
Reynolds number increases for flows around a cone, the flow asymmetry becomes strong and 
changes sides in the downstream direction. For the high-Reynolds number flows, the spatially 
asymmetric flow develops in a wavy manner. The mechanism of vortex shedding is qualita- 
tively similar to the temporal asymmetric flow of the locally conical flow solution, where the 
flow asymmetry develops in a periodic manner. As the angle of attack increases, the flow 
asymmetry becomes stronger and unsteady. The unsteady asymmetric flow case shows 
evidence of multiple small-scale vortices moving along the body and vortex shedding at each 
section. 

Adding a cylindrical afterbody to the conical forebody strengths the flow asymmetry in 
comparison with that of the isolated cone. Finally, a comparison of the computed surface- 
pressure coefficients with the experimental measurement for a cone-cylinder configuration is 
given. The results show that a slight flow asymmetry starts close to the nose region and 
spatially grows moving downstream. The shocks on the top of the primary vortices show 
strong asymmetry and change sides in the downstream direction. 

The next step for the research work is to study the control of three-dimensional asymmetric 
supersonic flow using a passive-control method in the form of side-strakes and/or an 
active-control method in the form of blowing or suction ports with various blowing rates and 
orientations of the ports on the body surface. 
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